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REMARKS 

Claims 233-247 and 466-479 are all the claims pending in this application. 

The pending claims are rejected based on 35 U.S.C. § 112, second paragraph, as being 

indefinite. 

Claims 234-235, 246-247, 466-467, and 478-479 are rejected under 35 U.S.C. § 103(a) as 
being unpatentable over ILIADIS et al. (IDS ref: Computers and Biomed. Res. (2000), vol. 33, 
pages 211-226). 

Claims 234-242, 244-247, 466-474, and 476-479 are rejecter under 35 U.S.C. § 103(a) as 
being unpatentable over ILIADIS et al. (IDS ref: Computers and Biomed. Res. (2000), vol. 33, 
pages 21 1-226) in view of PERIERA et al. (IDS ref: Frontiers Med. Biol. Eng. (Ie95), vol. 6 (4), 
pages 257-268). 

Claims 243 and 475 are rejected under 35 U.S.C. § 103(a) as being unpatentable over 
ILIADIS et al. (IDS ref: Computers and Biomed. Res. (2000), vol. 33, pages 21 1-226) in view of 
PERIERA et al. (IDS ref: Frontiers Med Biol. Eng. (1995), vol 6 (4) pages 257-268) as applied 
to claims 234-242, 244-247, 466-474, and 476-479, above, and further in view of LEMELSON 
(US 5,919,135). 

The Applicants traverse the rejections and request reconsideration. 

The Applicants have amended the claims to overcome each of the section 1 12 rejections. 

All the pending rejections are based on ILLIADAS, either standing alone or in 
combination with other references. As evidenced by the attached database report (Appendix A) 
from PubMed, the ILLIADAS reference, though having a "received" date of January 2000, has a 
publication date of June 2000. 
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The Applicants respectfully attach a declaration under 37 C.F.R. § 1.131, executed by Dr. 
Zvia Agur, swearing behind the ILLIADAS reference. The declaration and the attached Exhibits 
demonstrate a conception date of at least prior to June 2000 and continued diligence until October 
19, 2000, the filing date of the present application. 

Therefore, the Applicants respectfully submit that ILLIADAS does not qualify as prior 
art under any sections of 35 U.S.C. § 102. 

Since there are no other pending grounds for the rejection of the pending claims, they 
should be allowed. 

The Applicants respectfully note that the pending claims were first rejected in an Office 
Action dated May 8, 2001, based on the combined teachings of Barry (6,081,786), Fink 
(5,808,918) and Thalhammer-Reyero (5,930,154). After a response was filed by the Applicants 
on August 8, 2001, the Examiner maintained the rejections of these claims based on the same 
references in an Office Action dated October 25, 2001. The Applicants provided a further 
response on January 25, 2002. Based on this response, the Examiner appears to have withdrawn 
the rejections and issued the present rejection based on a new combination of references. 

In an unpublished, non precedential opinion, Administrative Judge Bahr asserted that 
"The application is remanded to the examiner for compliance with section 706.02 of the MPEP, 
which directs that prior art rejections should ordinarily be confined strictly to the best available 
art. . .Accordingly upon remand, the examiner should select the best prior art rejection for each 
claim and withdraw all cumulative rejections" (see generally, Appeal No. 2002-0965). 
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It is assumed that the Examiner has followed the requirements of MPEP §706 .02 and 
issued rejections based on the best possible references on both the occasions noted above; first, 
when the claims were rejected based on Barry/Fink^Thalhammer-Royero as well as, second, when 
they were rejected based on Illiadas in combination with other references. In this connection, the 
Applicants respectfully note that during the course of the prosecution while these claims were 
amended to overcome minor informalities, the subject matter of the claims has not been changed. 
Therefore, it is clear that the Examiner has already searched the subject matter recited in 
these claims twice and has found the best possible references to cite as prior art. 

Considering the above, it is respectfully requested that the pending claims be passed into 
allowance and the Applicants be awarded what is believed to be their rightful due. 



CONCLUSION 



In view of the above, reconsideration and allowance of this application are now believed 
to be in order, and such actions are hereby solicited. If any points remain in issue which the 
Exammer feels may be best resolved through a personal or telephone interview, the Examiner is 
kindly requested to contact the undersigned at the telephone number listed below. 
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The USPTO is direcKd and authonKd tQ charg£ au ^ eKept fa ^ 

P« a„ d the FM , , 0 ^ Accoum No ^ ^ ^ ^ ^ 

overpayments to said Deposit Account. 

Respectfully submitted, 




SUGHRUE MION, PLLC ■ 

Telephone: (202) 293-7060 Chidambaram Subramanian 

Facsimile: (202) 293-7860 Registration No. 43,355 

WASHINGTON OFFICE 





23373 

PATENT TRADEMARK OFFICE 

Date: March 18, 2003 
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MADE 



APPENDTY 

VERSION WITH M^ P KINGS TO sunw ni ^ ^ 

IN THE CF.ATMS . 

The claims are amended as follows: 

234. tended) A computer system for recommendl „ g „ optima] ^ ^ 
for treatt„ g cancer usi„ g drags , for an mdi ^ ^ ^ ^ ^ _ ^ 

said system further comprising: 

a cancer system model; 

a treatment p rot oc„, g e„e ra , or for „ g , ^ of ^ 
treating cancer using [chemotherapy] drugs; 

sa,d cancer system mode, based on parameKrs specfflc (q ^ ^ 
protocols based on the modified system model. 



240. tended, The system 0 f claim 238 where the system includes a se, of control 
functtons that are adapted to uniquely determine 



control functions comprise a g e of cells, state of 



an outcome of every sin g le step, wherein said 



a current population and associated environment. 



244. (Amended, The system of claim 243, wherein the system is adapted ,„ 
-corporate .pharmacoKinett, ph^^ ^ [pharmac0 „ Efcm ^ 
~ effects, cytotoxic effects, and other effects on eel, dilation „f micancer ^ 
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246. (Tw.ce Amended) The system of claim 234 wherein, said parameters specific ,o 
■he individual comprise parameters related ,„ , umor dynamics, patient specific dntg 
[pharmacokinetic, pharmacodynamic] e tmm^ mms ^^ a nd dynamics of 
dose-limiting host tissues. 



466. (Amended) A computer-implemented method for recommending an optima! 
treatment protocol for treating cancer using dn,gs [. induding chemotherapy,, for an individual 

said method comprising: 

creating a cancer system model; 

enumerating a plurality of treatment protocols for treating cancer using drugs; 
modifying the system model based on parameters specific to the individual; 
selecting an optima, treatment protocol from said pluraJity of treatment protocols based 
on the modified system model; and 

recommending said optimal treatment. 



472. (Amended) The method of Cairn 470 where a set of control functions uniquely 
state of a current population and associated environment. 



476. (Amended) The method of claim 475, wherein [pharmacokinetic and 
Pharmacodynamic,, E ta B a i c ia!! e fej ^^ cytotox.c effects, cytostatic effects 
and other effects on cell disintegration of anticancer dm. 



igs are incorporated into the model. 
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478. (Amended) The method of claim 466 wherein, said parameters specific to the 
individual comprise parameters related to tumor dynamics, patient specific drug 
[pharmacokinetic, pharmacodynamic,] ita^^ and rf 

dose-limiting host tissues. 
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Zvia AGUR, et al. 




PATENT APPLICATION 
IN THE UNITED STATES PATENT AND TRADEMARK OFFICE 

In re application of ^ 

Appln.No.: 09/691,053 ^ 1 

Group Art Unit: 1631 £L v 

Confirmation No.: 5359 lCv<s.\C15 

Examiner: Marjorie A. MORAN 

Filed: October 19, 2000 

For: sssrrrs™- 



DECLARATION UNT>FP ^rpp 



1.131 



Honorable Commissioner of 
Patents and Trademarks 
Washington, D.C. 2023 1 



Sir: 



I, Zvia Agur, hereby declare as follows: 

I am an inventor and an applicant of the invention entitled "System and Method for 
Op.imi.ed Drug Delivery and Progression of Diseased ^ ^ ^ ^ ^ ^ 
in U.S. Application No. 09/691,053 filed October 19, 2000. 

. graduated from the Hebrew University of Jerusalem and Universite Libre, Brussels 
Belgium in 1983, and obtained a Ph.D. Degree in Biomathematics. 

Presenfly ! am the Chairperson and CSO a, OPTTMATA, the owner of the above 

application. 

«or to June 2000, along with the other named inventors, had invented the invention as 
describe, and claimed in the above-identified apphcation, and pursued the present invention in 



•he ordinary course of business , , ncMng preparing above _ jdentjfied ^ 

the fiung of the above-identified appHcation, as evidenced by the fohowing: 

Prior to August 1, ,999, having earIier „ ^ ^ f ^ 

law firm of Eitan, P erl , Uaner and Cohe„-Ze« in Herzlia, fsrael, thereinafter "MC", to assist 

the applicants in preparing a patent application. 

On August 7, 1999, a meeting was held at Mc>s office ^ ^ 
he.d ,o convey the inventive subject matter to MC. Several related technical documents 
handed to MC, including the overall concept diagram as shown in Exhibit , . 



were 



were 



Subsequent ,„ the meedng on August 7, .999, several discussions were he.d 



over a 



period of several months while MC and the apphcants together prepared four manuscripts 
(Exhibits n-,4) disclosing the invention. Exhibits 2-8 taclu de tr a „s Cr ip,s of emails evidencing 

contmued work on the four manuscripts until June-July 2000. 

In July, 2000, we contacted the law firm of SUGHRUE MION, PLLC [hereinafter 

SUGHRUE] to provide further assistance and prepare patent applications based on the four 

manuscripts (Exhibits 11-14). 

Exhibit 8-9 shows emai, correspondence evidencing our discussions with SUGHRUE 
during July- August 2000. 

ft- ,„ August .8, 2003, technical materia, including the four manuscript. (Exhtbits , ,. 

14) were provided to SUGHRUE. 

We held a meeting on August 28, 2000, at the offices of SUGHRUE in Washington DC 
with Mess.. Bil, Mandir and Chid Subramanian to hold technical discussions on the subject ' 

matter of the invention. 




From the time of our meeting with SUGHRUE until the filing of the application on 
October 19, 2000, several drafts were exchanged between SUGHRUE and Optimata. 

The above-identified application was then filed on October 1 9, 2000. 

In view of the discussion above, and the attached exhibits, we had invented the claimed 
invention disclosed in the subject application prior to June 2000, and pursued the present 
invention in the ordinary course of business, including preparing the above-identified patent 
application, until the application was filed on October 19, 2000. 

I declare further that all statements made herein are of my own knowledge and are true 
and that all statements made on information and belief are believed to be true; and further that 
these statements were made with the knowledge that willful false statements and the like so 
made are punishable by fine or imprisonment, or both, under Section 1001 of Title 18 of the 
United States Code, and that such willful felse statements may jeopardize the validity of the 
application or any patent issuing thereon. 



Date: °? 




2. 



3. 



4. 



5. 



^ at a meeting on Aug 7 

E-mail from Dr.Agur to MC dated August 1999 con « . 

•Agur to MC dated March 16, 2000 
^en, search report aerated March 2J , 2000 

Mem ° 0faTO *8<ta tt dJune ,1,2000. 

pterin Hebrew dated July 21 2 0nn f „ 
Wween Dr.Agur aj un ' W0 fr0 ™ Dr.Agur to Mr j . 

g and MC regarding patent regi*^ C^^^ns 
E-mail dated July 31 20OOf ^ June -September 1999. 

yji, 2000 from Dr.Agur to MC. 

E-mail dated July 31 ?nnn f 

ui y Ji, zuuo from MrRiclr aku^ a- , 

9 p ., ClcAbboudl ofMCtoDrAeur 
E -^I to Bill Mandir with SUGHRUE d„ „ . 

10 p m ,. ^HRUE dated August 1,2000. 
• E- mai]s ln August STJrHprTr7 

11 T , . SUGHRUE and OPTIMATA 

19 _^ • yiCand Applicants. 

Document titled appendix A . r 
H r, fiINDK A. Computer Simulator of Human Th u 

13 - Document titled APPENDIX R * » Thr °mbopoi es i s . 

p^eunient titled APPENDIX C- n ,• . SF> 
(Ministration/efficacy) fofcSce" pa™^ W " ic *« 



7. 



8. 



St.. 



' ' • t » . t 



Bone Marrow Project- VBM 




r 



OPTIMA 



project 





VBM+OPTIMA 



project 





To: cohenm 

From; agur@ccsg.tau.ac.il (Zvia Agur) 
Subject: advice (Ely Lilly) 

Bcc: 

X- Attachments: 
Dear Mark, 



«ssk rp^^s^ r,r loowng — - - - 

for medical biomaths (IMBM) Thtll T I P portfo,io f <* the institute 

potential collaboration with Eli 2 y wni^T kT be ° f h ^ h ^ance " 
have to give the guy so me ini tt « m wiT 1 ? T 1 6nC,0se Wow. 
any suggestions about how to say to hJ^ZmZ I WOTder Whether ? ou »ave 
us enough information; (2) unlike whThf } CM d ° what he need « if he gives 
questions reliably is a sen s "ot^^T ** " Uw « i '* his 

operation research and computeTsdenV^ ^ h biol ° gV ' biomaths ' 

case I do not hear from you I w 1 answer h?J '1 m ° delS 3nd * is our * 
committment on my part y ^ ha ' ^ "^"J? « ^ *" ^ addl «g that any 
Zv ia Agur lo De a PP r oved by may attorneys. All the best, 



F±'.T" e ' 2 °/, ul 1999 13:16:51 -0500 



— , — ,u, ij.io.m -0500 
f rom: Craig M. Zwickl" <zwickl craic mffllT •« 
Modelling of Cancer/Hemato^ 
Ration: Eli Lilly and Comp ny ^^•tau.ac.il 
MIME- version: 1.0 



Dear Professor Agur, 



Ll " V and ^any in Ihdian^KAO' ,* ^ To ^ g ist at Eli 
and ecology at the Indiana Universk^ School nf 11* in P ha ™acology 

pnncmle investigator in a tema^Z^^l^T- M UUy > 1 am the 
effects of drug candidates on th^^J^"^ Which ch ^cterizes the 
colony-forming assays. Since^^S^ * vitro 

oncolyt, agents, I wor k f ai rly WSSJiX^ 

attributable to mathematics airarw^^J^ ^I^ 6 t0 US ( P robabl y 
that mathematical modeling of doTZl 7 M SCh ° ol) - Howev <* I am convinced 
gainmg an understanding i^^TjT be V ^ benefic -l for 
compound's therapeutic fndex thereby Zltn etc - ™uld affect a 

treatment. For this reason, I have ^ deveJoned , \ the effect 'veness of 

pnndples of mathematical rj^^^f h ™ to apply 
where to begin, as I have not bccn l^^J^^ bUt ^ know 
Your experience and interest in modeling ZT > V 8 ° d reSOUrces for self-teaching 
growth, and cancer therapy, n^^^ZT^**" ^ C ° ntro1 ' cance 
acquamted with the mathematical S^^^nSS^ * ^ me become 

quireo to model these systems. I am seeking 



*°™the opportunity t0 dis cus , his .n 4 tsarea 

U this is something which akn m ;„i,. • 

<c-®im,c„ m) . ihaCr y ?:z es ,'s ?™ - - «*< 

Regards, *° ^ X". 



Regards, 

Dr. Craig M. Zwickl 
Senior Toxicologist 

Lilly and Company 
PO Box 708 y 

Greenfield, IN 46140 

Date: Thu, 05 Aug 1999 10-08-55 Osnn 
From; Crais M 7»,^, " 050 ° 



SSS^r* Mathematical 



X-Lotus-FromDomain: LILLY 
Status: RO 



Dear Professor Agur, 



Thanks for your willinene« t„ w 

chemotherapy and hema ^Tj "»«■« » f cancer 

that I work on fhl r f f lnstructi °n you could offer Thet S P b,CmS and Wou 'd 
requires t T ^ ^ from ° f ^™ 

project. medl - ~ *■ t0 ^e-constramts l£ ^ 

At the moment, I can only discu« th* u, 

egress inth*' £ % ^^Str** 
compensation. ^'"^ ™» I mtgnt be abie „ ^atse"^^ 

Background: We have be™ eh ^ • 

"ZcZ T » teT*- Tte o" these 

P eaict plasma concentrations that would 



w^ptenS S^?"" I»«<™> ("I send you a C0DV of ,. 
approach assumes T"" 8 ^ k S&T'' 

gproach is the best for a „ Ita.^SSSSrS'- 1 ^ ^ 
tested, the ,„ v.tro inhibition curves for CmTr?/' ! raample - wi * the compound 

~d with cytoskeieta, struC It fe^T"* ^ *~ 

=^ESF-=^ore 
^^^^ 

hereby gtving a better therapeutic index H™ m ? ly ,ar8 « Proliferating^ 
e-nce te such . J2£%£ ^^J? 

S^tt^^SC ^ »» «« that have 
dependent coiony i„h lbition "= to <* very carefuliy at AUG vs concentaTn- 

Opportunities for Modelling if 

colonies at low kj^ZZZZT* ^ ^ 25% inhlbiti °" of CPU GM 

«P«ed ,„ lnfluence Ktjmales ob ,Le d C mo S„g? d ' UmOT '** "» 

I think these are very interesHnr, 

~ a= .,,.„„ „ 



compatible with your system? 

Regards, 
Craig 

P.S. I will be on vacation (much needed) 



next week. 



> 
> 
> 

> 

> 

> 

> 



2Via ' Patentmeeti ngl60300l.txt 

He was out last week and Monday j hoDP , 

Mark y ' 1 h0pe t0 s Peak with him this week 

Frnm~ 0 ^' 9inal Mes sage 

sen?.: £3? ZtffTffi'V.fr-W 

T °: Mark Cohen 0 9:14 AM 

Subject: 



> Dear Mark, 

> Lchty ed Very mucn our meet-inn i 



on the 
optimization 



Collaboration with Prnf m„ Ier apy. 
TJjjb very »uchtr y °„ f u Xf p - ^ value to us . 

°° k1 " 9 fon « a '' d t° hearing f rom you . 
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TOTAL FOR ^ FILES 

H26672 MATHEMAT? f 5A , „ n „ 
_ « «-0i, «'^« "S^JJ^ 0R DYNAMIC SIBW 

L7 26 FILE WPI NDEX tem Cell #) 

1 FILE JAPIO 

TOTAL FOR ALL FILES 
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- Original Message- 
Sen 0 ?: JJ¥y r g^?J?".or»J 

T°; Ricki Abooudi 00 2:45 PM 

sSbject'f^ C0h R e e n ! £arl zeev 
Rowing our phone conversation yesterdav , 

3 % expect the draft of the patent in thi , ^ ** ^ 1 
in JeSJSf"; ^ 1 -dose the fi le v ^ ^ JUly 

^via Agur 

Frnm^ Original Message 

^c^?£lOH"- T ™ C0 - IL> 

> > Dear zvia, 

> r^dy™ Sti11 Wrk ^ on the draft, and j will , 

^ > I was wonderinq w hpi-h»n you know w "en it is 

> > Regarding your w-ich +■ ■ 

Pagl l am str1 «ly the 



> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 



> 
> 



> > science/technology contact DeES^X 0700 -™ 

> > y> ^ UMLclt -t person at the present time. 

> > Sincerely, 

> > 

> > 

> > Ricki Abboudi 

> > 

> > 

> > 

> > 

> > Ricki Abboudi " 

> I SToff^' LATZER & COHEN-ZEDEK 

> > 2 Gav Yam Center 

> > 7 Shenkar st. 

> > Herzelia 

> > Tel: 972-9-970-9401 

> > Fax: 972-9-970-9001 

> > 

> > 
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Dear zvia, 



patentmeeting310700.txt 



i am sti I I workinq on the draft t . -n n 

I was wondering wfiether you ?oJid Send 12 LT y< ?Y k 2 0w when ^ is ready. 

appendix B (g-csf model) the wav In?, con* y e " mai1 J he text of your 

and c (cancer). U ' the way you sent m e appendix A (thrombocytopenia 

Regarding your wish to dkmcc 

protect Wur patents their co^fS a l S SfS^ff nCernin9 the abilit V to 
legal/commercial questions to Mart r«L 9niflcance ' etc -> please direct a 

Sincerely, 



Ricki Abboudi 



Rlcki Abboudi 

M™6ffice RL ' LATZER & C0HEN -^DEK 

2 Gav Yam Center 
7 Shenkar St. 
Herzelia 

Tel: 972-9-970-9401 
Fax: 972-9-970-9001 
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Subramanian. Chid 



Sen?"' Moshe Vardi r m °shv@imbm.org] 

To: ' IZ^T* 20 ' 2000 8:5 8 AM 

Cc . Mandir, William H. 

Subject- CSubramanian ©sughrue.com 

ouojeci. Re . Qptimata meeting 



Dear Bill, 



2K KSrSt* W6 haVe talked and 
therefor you can feel free Tn^^Z^JTT' ^ 
everything on the 28th (9:30 will be comEntem > J? u 3 " dlSCUSS 
say that you have a green ligh o nroS 29th " 1 would 

to bring us to a poin 'where ne^Sc ml ^ Ur JOb ' and 1 ex P ect this 
"getting to know each other" bu S W ° U ' d be more than a 

question of how you and Mr SuhSS S ln my op,nion is mo s% a 
here this week toVnswer o ' ^M^iSSS ^ the 1 am 

Last but not least, it will be a o llSnr^n h V te ° hnical lssue S- 
would like very much to ?hank Sf w tL ^ d,nner with ^ and ' 
is less dangerous'n te?ms of uZak^r ^SJ 0 "- Mondav the 2 ™ 
Sunday, but I cannot giarentee fol S ' me 0 " be landin 9 on 

I m looking forward to meeting you and Mr. SubraTaS ?" 

Moshe 

- — Original Message — -- 

M f d /[' William H - <wmandir@sughrue com> 
To. Moshe Vardi' <moshv@imbm.org> ™ 

Spnt SL ama A ian ' Chid < c subramanian@sughrue com> 
Sent. Friday, August 1 8, 2000 9:39 PM U9nrue - com> 

Subject: RE: Optimata meeting 



> Shalom Moshe: 
> 



> on Monday aSoon Xtl?: ,a J that vou "ave a meeting 

> 2100 Pennsylvania Avenue Nw Jhf^e ?„ Ur ° fflCes around 9:30 am at 

> Pennsylvania ^ a TSe SSSl^SK £e' "2? 2< 

jw»y close to me Foggy Bottom Metro stop. There w» be a directory in 

>' S^'^^JS; £f. -V, our drees code Is business 
of s ' Al ,hls P 01 " 1 ' 1 ar " "nable lo gauge the amount 

> meeu am ava,,ab,e a„ o^nTa^^ Ite <° 

>xrsa»^^ 

> TechrCaXs'Se^hKSr *""! 0,6 "» - 

> Science. Chid alsota fa S JS^ ?wT"?' en 9 |nee *9 Computer 

> puestions-nis *. has /ffiTSK^fflrSS'^ 
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> necessary conversion software. Therefore, pleiTsltie Regain 
» a WORD formal directly to Chid at CSubramanian@Su 9 hrue.corn. 

: sees '&sss t'X'zfzTr Mon r 

tec n rS° a <a " d WednesIla V) * *> left open for furlher 

» SnToSoSS, S™"< P-P—V comments for our 



> Patent Strategy 
> 



> of VC money and was still refining the invention, so the option of filing 

> provisional application, which is somewhat less expensive in comparison to 
^regular application, especially if new subject matter will be added, is 

> appealing and made sense. On the other hanH it 

> Fortune 500 company, and t^lZ^^iS^S^S!! WaS f 

> perhaps a provisional application would not make m ?rh con ° d ® velo P ment . 

> mentioned that time is crrtirai L ihi* „»?l, uch sense - Your e - m aN 

> co^pTny 6rm ' ' 3m C ° nfident We Wi " devel °P the strategy for your 



> 

> Communication 
> 



t aSked about wa y s of communicating. I can tell vou how l 

> ^! ust ° ur scnedu| e knowing that Israel's work week is 
»^«a»^^ 

> will visit Israel next month 

> clients-in comparison, communicating with Israeli clients is a breeze! 
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> World-Wide IP Coverage 

> 

> We can advise you as to how to protect your invention in the U S 

> and around the world. We have long-term relationships with foreiqn 

> associates all over the world, and we have our own office in Tokyo Japan 

> We cannot guarantee coverage in the U.S. or anywhere else in the world--! 

> wish we could, but none can. What we can do is put you in the best 
position 

> to obtain the maximum possible patent coverage. Every patent in the U S 

> and abroad requires examination by a patent examiner. For that reason ' no 

> one can possibly guarantee coverage. However, we have vast experience in 

> dealing with the U.S. Patent Office, as well as patent offices SSZx 

> the world, and based on that experience we can advise you so your 
invention ' 1 

> is in the best possible condition for allowance. 



> 

> Costs 
> 



> Based on our preliminary examination, it appears your disclosure 
prepare' 3 C ° nta ' ns Several inventions (perhaps 3 or 4). Assuming we 

> one application containing all inventions, we would estimate the cost to 

> ?f, a nS?" a PP |ication in tne U.S. to be approximately $8,000 to 

> $1 1,000, including all government filing and other fees. A provisional 

> application would be significantly less, about $2,000-$2,500 

> We look forward to meeting you and Prof. Agur on August 28 

> Finally, we would like to invite you and Prof. Agur to dinner either on 

lee 27 ° r M ° nday " AuQUSt 28 - Please let me know if vou are 

> and if so, of any dietary restrictions. I know of a terrific Kosher 
French 

> restaurant but it is difficult to get reservations, so let me know if you 

> are available. I will be out of the office until August 25 but will be 

> checking my e-mails. Regards, Bill 
> 

> William H. Mandir, Esq. 

> Sughrue, Mion, Zinn, Macpeak & Seas, PLLC 
>Tel: 202-293-7060 

> Fax: 202-293-7860 

> 

> > — Original Message — 

> > From: Moshe Vardi [SMTP:moshv@imbm org] 

> > Sent: Friday, August 18, 2000 1 :26 AM 

> > To: William H. Mandir 
Subject: Optimata meeting 



> > 

> > 

> > Dear Bill, 

> > 

> > 



Myself and Prof. Agur are in Washington on Aug. 28, 29th. We have 
scheduled a meeting for the afternoon of the 28th in the city, and would 

> > very much like to meet with you and your team prior to that. I would 

> > consider staying in Washington any needed additional time, to enhance 

> > patent filing procedure. I urge you to consider the opportunity to move 

> > forward the patent drafting, taking full advantage of my stay in 

> k > Washington, and the option to delay my flight back. I would very much 

mind S ° me W ° rk W ° Uld be d ° ne by that time ' and through a day or two of 
>> storming we could bypass all possible obstacles. Nevertheless, I cannot 

> > so without your advise. 

> > Please reply ASAP. Note that Prof. Agur is leaving Israel on the 22nd, 

3 
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FIELD OF THE INVENTION 

an, cor ' nVen,i0n relateS '° 3 SySt6m - * 

T °' '~ - *~ « - certain s,a 9 e s , n dru 

development and trials. 9 

BACKGROUND OF THE INVENTION 

T ,ha ' are adminis,ered ,o combat dfeease - - - «. 

~ th ee ff e«o ftheprotoeo , onbolhabnormaiandnoniia:cei(s 

As a JT*' m0de ' S bi ° ,09fca ' ~ ^ ^ - « 
a broad cate 90fy, these models are... 

More specifically, can have models of cell l,w , 

a or ce " "nes, tumor growth etc 

^^^^^^^^ 

;. aPPr0aCh ' ** ^ — -~ —on, o anno( be used t0 

::r of,r ~ inrea - — — — 

«*-» o, more complex and detailed scenarios is needed. 



SUMMARY OF THE INVENTION 

completed when claims are finalized. 



BRIEF DESCRIPTION OF THE DRAWINGS 

The present invention will be understood anrt * 
from iho f it • naerstood and appreciated more fully 

from the following detailed description ^ • 

description taken in conjunction with the 

drawings in which: 



appended 



understanding Fig. 1 ■ 



Fig. 3 is a schematic illustration 
with one embodiment of the 



10 



different format; 



of a biological model, in accordance 
present invention; 
Hg. 4 is a char, iteraaon of the bioiogica, mode| of Rg 3 . 

5 is a graphical iteration of ,„ e chart of Rg . 4 . 
R3. 6 is a char, iteration of the biologicai mode, of Fig. 3 in . 
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of Fig. 3; 



R ' 9 - 7 iS 3 9raphical Nation of the chart of Fig. 6; 
Rgs. 8A and 8B are graphic illustrations of 



the output of the model 



experimental data as 



and 8B; 



model, in accordance 



Figs. 9A and 9B are graphical illustrations of 
compared to the output shown in Figs. 8A 

Fig- 10 is a schematic illustration of a biological 
-th a further embodiment of the present invention; 

model shown in Fig. 10; 

Rg, , 2A and 12B are graphfca| taratas ^ ^ 

*se S o, G -CS Fo n,h sneutrophil , e , ao:ordingtothemodeiofRg ^ ; 



Fig. 13 is a schematic illustration of a biological model, in accordance with , 
further embodiment of the present invention; 



INVENTION 



DETAILED DESCRIPTION OF THE PRESENT 

A system and method have been developed for identifying optima, 
-atment protocois for se,ected parameters. Un „ k e prior art option 
-hemes, the present invention uses b*^,, The ^ 

5 heUnS,i0S enaWeS ,0 " ^JSSSLSSSUm * conpte* (h enoe 
-istio, mathematical ^.^ ^ ^ ^ 

case ("optimal generic treatment"), as well as at th* ^ . , 

i, wen as at the level of an individual 

P^ent Reference is now made to Fig. 2 , whioh is a flow ^ ^ 
-o general showing , he optimization method. The Figure genera,,, depicts *e 

Determination of op„ma, protoco, is cased on specKic parameters input by a 
The user may he a phys,ia n , a drug deve,oper, a scentist, or anyone 
« else who may need to define a, reatmen, protoco, for a drug. The specific 
Parameters may inciude several categories, such as ind^dua, patien, 
characte„s, ics anoVor med,ca, history, needs o, a specific user (research 
erccacy, treatment, etc, and other particulars (such as maximum length 0 , 
treatment, confidence level, etc.). 

The firs, steps include generation o, several types of models. A se, of 
-dels for a„ the relevant hiologica, processes, is created, in addition a 
-de, o, treatment effects on each o, these processes is created The 
combination o, these mode,s proves a de,a„ed mathematical mode, of the 
overa,, bicclinica, scenario in a genera, or specific sense together with the 
specific effects of a part,cu,ar treatment. Once the comprehensive mode, is 



~e d the earache parameters (eNher „,„ ^ ^ 
P-n, specific e.g., - age , gender , we . h , ^ c|injca| . dicatfons) 

fomented in , ,„ tnjs way , ^ ^ ^ ^ ^ ^ 

definition of a p rot oco, space. To do this, possib,e permutations of certain 

5 Parametera " ^ *~* *** «-a,s, drug quantlty , 6tc are 
considered. Thus, . number (can be ^ iarge) ^ possiWe ^ 

P-ocois is generate* The amount of possMtfes ^ ^ ^ ^ 

ranges of parameters considered. At this point, , ne fitness function is 
constructed by mathematically ^ possibie ^ ^ 

- may be fenced by the treatment. These may include survival, tumor or 
Pathogen ioad, cytotoxicity, side effects, pain etc. The user can pu, certain 
specific parameters, the «ness function so as to a dj us, this function to 
«her specific goaisBased on the seiected parameters, the fitness function is 
applied, and calculates a fitness score for each and every possibility in , he 

> Protoco, space. Rna „ y , the „„„ step fe ^ ^ ^ ^ ^ 

heuristics or by analytlcal methods , ^ ^ ^ ^ ^ ^ 
protocol from all the scored possibilities. 

«-» * — once the parameters are entered may be negiigibiy short or up to 
hours, depending on the length o, me simulated treatment period and the 

Power of the specific search heuristics and the computational tools, making 

this a very feasible tool. 

Construction of detailed mathematical models for biologica, processes 
and treatments „i„ be shown by way of example. A mode, o, piatlets 



production and thrombocytopeniawith treatment by TPO, a mode, o, neutrophil 
production and neutropenia and treatment by G-CSF, and a modei of cancer 
growth and chemotherapy are described. An exampte for specie opUmization 
»y linear programming, i s imp teme n.ed in the neutrophii mode,. A genera, 
5 heuristic optimization method is described as well. 

EXAMPLE 1 : Thrombocytopen la and Treatment by Thrombopoietin (TPO) 

Thrombocytopenia is a common hazardous b,ood condition, which may 
appear in different Cinioa, situations, inCuding cancer chemotherapy. 
Recently, a thrombopoiesis-controlling cytokine, thrombopoietin (TPO, was 
' *°>a«ed and „ s human recombinant analog became available A 
mathematical mode, has been developed to imitate dynamics of a 
thrombotic lineage in the bone marrow, of p,a,e,e, counts in the periphery, 
and effects of TPO administration on them [on what? Platelets?]. 

TPO is a cytokine, glycoprotein of about 350 amino acids that 

resembles erythropoiesis-stimulating hormone, erythropoietin. ,,s synthetic 

ana,ogs, recombinant human thrombopoietin (rHuTPO, and recombinant 

human megakaryocyte growth and development factor (rHuMGDF) are 

avai,ab,e as we,, and are undergoing Cinica, trials. These compounds have 

been shown to have the same biological activity as TPO has, so the term 

TPO wi„ be used without distinguishing between its different forms and 
analogs. 

TPO is a primary growth factor of the thrombotic cell line both in vivo 
and „ vi.ro. Aside from this, TPO may be potent in s«imu,a.ion and co- 



—on of other hemopoietic , ineages ^ granu|opofe(jc ^ 
erythropoietic). 

A. Background of Thrombopoiesis 

Like a„ other hemopoietic lines, the thrombopoietic line originates from 
poorly differentiated, multipotentia, cells, capable o, some division and se„- 
reconstitutio. Other b one marrow cel, iines, such as potential 
hemopoietic stem celis and common myeloid progenitor ce„s (CFU-GEMM) 
have these characteristics as well. 

' Starting ou, as undifferentiated stem cells, the cells gradually become 

-re and more differentiated and committed to a specific line. Once they 
commit to the thrombopoiesis iine, they proliferate extensive, within certain 
oompartments. One such compartment is tha, of colony-forming units - 
megakaryocytes (CFU-Meg,. other examples of compartments burst-forming 
™<s - megakaryocyte (BFU-Meg) and megakaryoblas, compartments, which 
have similar properties as one another. 

The committed thrombotic cells, which are megabyte 
Posers, go through several stages of maturation when proliferation is no 
longer possible. However, megakaryocyte maturation fe somewhat ^ 
from maturation within other hemopoietic iines. In megakaryocytes, cel, nuclei 
undergo mitosis in paralie, with cytoplasmic maturation. However, although 

incomplete mitosis is called endomitosis. Conseouentiy, the cei, becomes 



chromes are called either megakaryoblasts or megakaryocytes. 

Usua,,y, megakaryocytes do not star, , 0 r e,ease piatelers un,„ the ,6N 
state. Then they begin to create demarcation membranes tha, enveiop 
' cytopiasm fragments generating plate(ets [wha , does thfe ^ ^ ^ 

cell volume and an ability to release niaMat* • 

y release platelets all increase proportionally to the 

cell ploidy. 

B. Mathematical Model 

10 Re,6renCe iS "° W ™* ,0 * 3 ' *ch * a detaiied lustration of a 

m ° de ' PrediC,in9 ,hr0mb ~ A * * B, 3, , h e thrombotic 

"neage i s divided into eight compartments. The firs, compartment, caiied 
Stem ce„s (SC) and iabeied 30, refers to a„ bone marrow hemopoietic 
Progenitors that have an abiiity , 0 difference into more than one „ne (e g 
* Plunpotentia, stem ceiis, CFU-GEMM, and others, Ce„s of SC compartmenl 
30 proiiferate, mature, and subseouentfy differentiate into melanocytes or 
Cher precursors, or they may give rise to new stem ce„s. Afthough the term 
"new" is no, bioiogicaiiy accurate, , se„es as an accepts assumption for 
the purposes of the model. 

Ceil dea,h .hrough apop,osis may have a significan, effect on cel, 
number wi,hin proliferating compartments. However, the effect of apoptosis is 
combined with the effect o, cei, probation into a total amplication of cel, 
number in a given compartment. Thus, for exampte, the amplication rate o, 
cells ,n SC compartment 30 is normaiiv , , 029 „ew celis per hour An 



assumption is made that no apoptosfe ^ ^ 
megakaryocyte counts, due ,o ,ac k o, evidence ,o the contrary. 

Biologically, rates of proration and maturation, the ability to 
— tute, and other charac,e„ s ,ics differ between particular cei, types 
> within a primitive progenitor popuiation. However, in this model thsre „ no 
distinction between then, a„ progenitor cells are considered to be one 
population with common properties. 

Based on some works (which ones?, , ha , have shown , ha , 
of stem cel, differentiation into one or another hemopoietic ,i neag e are 

assumed about the stem cel, se, f ,enewai. Thus, after the cells spend a 

defined transit time, for example 24 hours in c?r . ^~\^ 

nours, in SC compartmen(30) a certain 

instant fraction of the celis ,0.5,81, return to their >ung state', i.e star, 
' 'heir passage through SC compartment 30 again, ^^1^ 1 
Another constant fraction (0.34 2 , for example) of otH^^T^ 
compartment calied Co,ony-F„rmi„ g Units (CFU-Meg,, labeied 40 „ ,s 

P~ tha, remaining stem ceiis differentiate into other hematopoietic 

lineages. 

CFU-Meg refe. to a„ cells tha, are already committed to the 
mega^ocyte „ne bu, are still capable of proliferation. Cells of CFU-Meo 
compartment 40, like those of SC compartment 30, spend some time 
(normaily 60 hours, multiplying a, an amplification rate of about t.0255 (per 



wha, ? , and maturing before losing their probative abiiities and passing on 
to the next compartment 50, called megakaryoblasts (MKB). 

MKB compartment 50 includes a„ , he c6 „ s tha , have lost the ability to 
proliferate, but are no. ye, sufficiently mature to release platelets. For the 

* purposes o, the mode,, the assumption is made tha, megakaryocytes do no, 
star, ,o release platelets until ,hey reach ,he ,6N-p,oidy phase. Hence MKB 
refers to 2M, 4N and 8N cells o, megakaryocyte lineage ,ha, canno, divide, a, 
all stages of cytoplasmic maturity. After these cells spend the designated 
transit time (normally 132 hours, in MKB compartment 50, they move to the 

» next compartment 60, which is an MK16 bone marrow compartment. 

The cells of MK,6 compartment 60 are megakaryocytes of ,6N-p,oidy 
class that release platelets until they exhaust their capacity, and then are 
disintegrated. Cel, voiume has a linear relationship with megakaryocyte 
plo,dy. Hence, we assume that al, ^-megakaryocytes have the same 
volume and, thus, the same platelet-releasing capacity. Therefore a„ platelet- 

releasing 16N-megakaryocytes are in transit fnr iho „ 

/ j a aic in transit for the same amount of time 

(120 hours) until they are exhausted and disintegrated. 

However, some ^-megakaryocytes do no, participate in platele. 
release, bu, rather continue with another endomNosis over a 48-hour time 
period, and become 32N-megakatyocytes. These constitute a new and 
«inc, MK32 compartment 70. Thus, after 48 hours in MK16 compartment 
60, a certain fraction of the cells (normally 0.26) leave MK16 compartment 60 
and go on to MK32 compartment 70. 



32N.megakaryocy.es re|ease p|ate(ets ^ ^ ^ ra(e o( ^ 
release is constant tor eve,y compartment and proportional to the pioioy state 
of megabytes in it. Thus, every , SN-megakaryocyte releases ,0 2 
Platen per hour and eve* 32N.megaka oocyte releases 20.4 p,a,e te ,s per 
' nour. However, 32N-me 9a ^ocytes are no, exhausted more ouicKly than 
^-megakaryocytes, since ,hey have 2 ,imes grea.er volume and piateiet- 

releasmg capacity. Consequently all platelet • 

*' a " P |ate let-releasing megakaryocyte 
compartments have the same transit time. 

Once again, some traction o, cells (normally 0.0,6, are no, engaged in 
-o Plateie, tormation, and con,inue to the 64N-stage. Additional endomitosis in 
MK32 compartment 70 .akes the same amount ot «me ,48 hours, as in MK ,6 
compartment 60. The 64N-megakaryocy.es continue the process in a new 
MK64 compartment 80, and 0.4 o, them become ,28N-ce„s in ye, another 
MK,28 compartmen, 90. Additional endomitosis in MK64 compartment 80 
15 lakes the same amount of time (48 hours, as before. 

Norma,,, 64N- and , 28N-megakaryocy.es are no, found in humans 
conditions. Therefore, we inciuded these compartments in our mode,, bu, in 

ovocytes ot greater pioidy classes have no, been encoun,ered in 

humans. 



Rnally, there is a platele, compartment ,00. This is no, a bone marrow 
compartmen,, bu, ralh er the piateie, poo, In the peripheral b,ood f(spIe e„ 
sequestration^,,. P , atelets refeased from ^ ^ ^ 



64N-, and ,28N-p,„idy classes enter p, a ,e,e t compartment 100 . Platele, 
elimination from ,he b,ood is presumed ,o occur after a defined .ransi, time in 
c,rcu,a«ion by a mechanism of death from senescent [wha « is To 
simplify the mode, we ignore the destruction of piateiets by their use [■>] 
s Transi, time through platelet compartment ,00 is 240 hours but can vary, as 
will be discussed below. 

II. MODEL OF TRfat^fmt rriTm 
A. Background of TPO 

The major sites o, TPO production are the liver and kidney. TPO is also 
» produced in the spleen and bone marrow, bu, the production rate in these 

and kidney TPO producion are cons,an, under thrombocytopenia and 
.hrombocytosis o, varying degrees of severity. TPO production in the spleen 
and bone marrow is inversely related to the megaka^ooyte mass, bu, ,he 
actual contribution is negligible with regard to total TPO production. 

Another mechanism o, TP0 concentration regulation is receptor- 
mediated TPO uptake, since TPO-receptors on the platele, and 
megakaryocyte surfaces are ,he main TPO-ciearance mechanism. Thus 
TPO concentration is inversely related to the total plateie, and megakaryocyte 



mass. 



The effects of TPO on the thrombopoietic line may be divided into three 
types: (i) stimulation of proliferation of megaka oocyte progenitors tha, have 



an abi,„v to pro(iferate; ^ ^ ^ ^ 

P~; „ induction of addi(iona| endomi(osis ^ a(re ^ ^ 

megabytes, which ,eads to an increase in tha moda| 

nlnirlv/ 7 



B. Mathematical Model of TPO Effects 
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TPO concentration effects on the thrombopoiesis line is now considered 
As discussed above, , hre e thi „ gs depend 0 „ TpQ " 
amplification rate ,am p) , tne rate of ce| , ^ ^ 

<""• *** a given compartment (transit,, and tne fractfon „ f 

Oocytes ot given pioidy tna, undergo additJona, endomitosis and pass 
on to the next ploidy class. 
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1 ■ TPO nnnro n t mt j nn 

Recombinant human full-lenath TPn ^ « * 

lengtn TPO and its truncated form rHuMGDF 

are fully active biologically. Wore, in our mode, we add exogenous,, 
administered recombinant protein to endogenous,, produced TPO in order to 
calculate actual TPO concentration. 

As mentioned above, the rate ot TPO production in the main TPO 

Produce ' *" " "** ^ — « under thrombocytopenia or 

thrombocytosis. The level of TPO mRNA in sites like th» h„ 

sues UKe the bone marrow and 

- no, significant, different from the TPO ,eve, in perfphera, biood. Therefore ' 

- assumption is made tha, the bone marrow and spleen contributions to the 
total TPO concentration are insignificant. Endogenous,, produced TPO is 



assumed to have a constant rate of production of 7 pg/m ^ ur ^ (his 
number can change. 



The main mechanism , hat controls TPO concentration in the biood is 
5 re ° eP,0r -~ d ™ Both megakaryocyte and piateie, mass 

contribute to the tota, receptor number (normaiiy 8 . 375 x ,0°, and 
oonseouentiy, to the rate o, TPO Cearance. We assume that every piatele, 
bears 220 TPO recepto. on what? , Each ^ ^ 

that can be reieased dunng the rem a ining rrfe-span, times 220 receptors per 

potential platelet. We assume «, . 

we assume also that every receptor molecule removes 

4.776 x 10' 10 pg/ml/hour of TPO from circulation w„ 

^ irom circulation. However, this number can 

change. 

Another mechanism of TPO remove, from the blood is non-specitic 
• TPO-receptor-independen, Cearance. This mechanism is rather indicant 
-n <he norma, s, a ,e, when rece P ,or- m edi a ,ed TPO binding, endocytosis end 
degradation remove most o, the TPO. However, when the amount o, TPO 
nses signiticantiy ab0 ve the a bi,ity o, the receptor poo, to remove i, non- 
specific Cearance becomes important. ,n our mode, this type 0 , TPO 
clearance is exponent, i.e. every hour some traction (0.03) of . current 
amount of TPO is removed from circulation. 

,™.ia, maximum TPO b,ood concentration to the administered intravenous (IV, 
dose (the relation coefficient is 0.0167, based on literature data 



The state when TPO completely disappears from the Wood is very 
unlikely, , he iower limi, o, possibie TPO eonoen,r at ,on is res.ric.ed .o a certain 
minimum, in this case 0.01 pg/ml. 

Thus, ,he formula ,ha, calcula.es TPO concen.ra.icn hourly is given in 
5 Equation 1 as follows: 



C M - maxftq + Endo+ 0.01 67- Exo- Nc-C t -R c - Rp),0.0l) 



d) 



where C, is TPO concen,ra,ion a. .he curren. hour („, Cm is TP0 
■0 concen.ra.ion a, ,he nex. hour, a* „ endogenous* produced TPO per 
hcur, a, is .he ,o tel amoun, o, TPO administered imravenously during the 
current hour, * is , he coefficient o, non-specific TPO ciearance, * is , he 
coefficient of recep.or-media.ed TPO clearance, i.e. amoun, o, TPO .ha. each 
receptor removes per hour, and Rp is the total TPO-receptor pool. 

In steady-state TPO concentration is cons,an, and equals ,00 pg/ml. 
2. TPO efforts 0 n am plifirat.^n mtr 

In our model, there are only <wo compartments, SC compartment 30 and 
CFU-Meg compartment 40, whose cells are capable of dividing. These 
compartments differ signtficantiv ,rom each Cher, so ,hey wil, be discussed 
- separately. Cells of o,her compares do no. proliferate, so .heir 
amplification rate equals 1 under all circumstances. 

a. SC compartment 30: 

Cells of SC compartment 30 are relatively insensitive to TPO as 
compared ,o committed megakaiyocytic cells. In our model this is considered 



as a thresh0W ot 2500 pg/m , of Tpo m ^ (hjs ^ ^ 

7 ^ ~" - «* — * — ,lo, As :ong as 

7" be ' 0W * - «* ~ W * an ,„ Mnsic 

T P 0-,nde P enden t mechanjsm tha , keeps sfee ^ ^ ^ ^ 



5 constant. 



Thus, below the threshold <;r ^ x . 

"resnold, SC ampttaanon rate (a mpsc) is delermined 

Ending on ths CUTOn( ^ Qf ceNs ^ ^ 3Q ^ 

dependence e quati on gi ves a sigmoidal funotion: ^ ^ ^ ; 
"0 amplfation , the ^ numbef ^ me ^ ^ (o ^ ^ 

when the ce,l number approaches inMy „ „ ^ ^ ^ 

<=- -ber is norma, [d e«ne „o m a„ ampsc equals one ^ , ^ 
— Pos Sib ,e value. S,em cel, arnpfca „ on rate „ ^ 
as follows: 
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am Psc (quant sc )= (a mp 

v r SC max 



-l] q uantnorm sr s 

3 • 5c 5 + quantnorm , n 



(2) 



""ere ampsc ^ is the maxima| ^ ^ ^ ^ ^ 
~en, 30, „ OTsc „ , he noma| ^ Qf ce|(s ^ gc 

~* - - ~ «s the arnpWication 

~, y due t0 sllen, changes of ^ ^ . ^ ^ ^ ^ ^ 

be. S genera equals 0.2, but this number c an be oh a nged. 



<n a range of very ,ow ce„ numbers (100 times ,ower then norma, or 
'ess), the mode, must define the stem ce„ parameters for the quickest 

recovery of SC compartment 30. Thus amn in ^ 

nus, ampsc in these conditions achieves 

the maxima, value (amp SCmax ). 
> On the other hand, the rise of TP0 concemration ^ (hreshow 

should occur in severe p lat e,e t and/or megakaryocyte deficiency, or when 
TPO has been administered exogenous,, ,n these situations the hemopoietic 
system mobiiizes addiiona, resources for recover o, piateiets and platelet- 
"'easing compartments. Hence, we assume ,ha, in an effective 
- concentration, TPO increases the rate of eel, ampliation in SC compartment 
30. This is shown in Equation 3, as follows: 



awPsciwantsoC^ 



( 



( am Psc^ -1) Wnor^ 



3 ■ quanta + 



+ 1 



sc + quantnorm 



l sc 



+ ( am Ps Cnm -l) F*(c'(C f -threslSt) 



(3) 
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where the firs, of two operands is the amplication calculated based on a 
TPCindependen, mechanism, and the second is the TPO-reiated addition to 
the amplification rate. 

G- is a function tha, transforms concentration (C, into a form tha, is 
easy to use in the calculation of both ampsc and tens/fec as will be dfecuss6d 
futher [where?]. G" is defined by Equation 4 as follows: 



G*(c) = ln(C + l). 0.23 + 0. 



03949 



(4) 



P- in Equation 3 is a function tha, determines the steepness, denoted * of 
-he a^versus-C, curve. F", as shown in Equation 5, is a recurrent funefcon 



•ha. adds 1 ,„ its firet argumen , and app(ies a iogari(hm (o thfe 

again. This operation is repeated St times. 

F '(*)=... log (log (log (log (* + 1 )+ 1 )+ 1 )+ ] ).. (s, 
Although S, appears in severe, equations, its value is specific ,or eveiy case 

•ranscends the threshold, the TPO-related addition is calculated relative to the 
difference between an actual TPO concentration and the threshold. 

b. CFU-Meg Compartment 40: 

15 

in contrast to the oelis of SC compartment 30, oelis of CFU-Meg 
compartment 40 are very sensitive to TPO and respond to the absolute TPO 
concentration, rather than its ieve, above a certain threshold. Biologically 
these oe„s have no TPO-independen. proliferation mechanism, and thus 

20 cease to proliferate when deDrived of Ton u 

aepnved of TPO. However, in our model, CFU-Meg 

cells continue to proliferate without TPO a, almost the same rate as with TPO 
With norma, TPO concentrations, ampm equals one eighfh of ife ^ 
Poss,ble value ( amfW The equa „ on ^ ^ ^ ^ ^ 

amplification rate of CFU-Meg cells < ampcFu) t0 TPO concentration <Q is- 



amPcm ( C ')=( am Pc mw -l)-F(G(c,\Sr) + i 

(6) 



Function G resets function G* (Equation 4) qua|jtatjve(y ^ ^ 
quantitatively, having the following form: 
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G(C)= 2.56 x 10 - 6 • (ln( C + 1)+ 0.3949 ) 7 + 0.3 



(7) 



(a) 



(b) 
(c) 



(d) 



of 7, multiplying by 2 , 6x10 -a and ^ ^ ^ ^ ^ ^ 
io following features: 

in the region near ,he norma. TPO concentrations the function 
rises with an increasing rate from about 0.3 up to about 0.8; 
in the normal TPO concentration (100 pg/m|, i, equals 0.5; ' 
in high TPO concentrations, which occur in cases of exogenous 
TPO administration, the function rises with a decreasing rele- 
ase (iii,, in high TPO concentrations tha, may arise in a 
logical situation (i.e. severe thrombocytopenia or TPO 
administration in logica, doses,, the rate of function growth does 
no. drop down to zero and the function does no. reach a 
maximum, bu. ra.her the rate gradually stabilizes a, a nearly 
constant value and the function continues to hse almost linearly 

Equation 5. The difference is tha, while F* returns the number sole, on the 

0.5, no matter what Sfis. This is achieved by multiplying the returned number 
by a certain coefficient, which is specific for every value o, St Thus, F is 
defined by Equation 8 as follows: 



F (*)- * • (... log (log (log (log (jc + 1) + i) + i) + !)..) (8) 
where kis the S/-specific coefficient. 
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As need earlier, a„ platelet-releasing megakaryocyte compartments 
have the same trans,, time. Since neither the r6lationship of megakaryocyte 
- volume (and thus, its p,a,e,e t releasing capac ity) to megakaryocyte pioidy nor 

compartments are constant, for example ,» hours. P la ,e,e,s in circulation 

also spend constant period of time ten ion h«, n u- . 

u or ume (eg. 120 hours) which is not affected by 

TPO concentration. 

in contrast, SC, CFU-Meg, and MKB compartments have changeable 
transit times. „ the trans* time were continuousiy strictly related to the TPO 
concentration, then during an abrupt concentration change (in exogenous 
administration, tor exampie, there could be either massive cel, exit from the 
otven compartment or absolute absence ot cel, exit tor the period ot time tha, 

station seems biological, unlikely, in lhe presem mode , ^ ^ 

changes gradually, approaching the value that it should be equal to based on 
TPO concentration step by step. 

Regarding the amplication rate, transit time in SC compartment 30 
25 differs from the other compartments. 

a. SC compartment 30: 



Cells of this concern respond ,o TPO only when its 

«*~ by an intrinsic mechanism d epen d en, on the ourrent C9|| number 
onl, The value tha, the transit time shou ld approach Is determine. based on 
•he cell number, according to the following equalion [=?]: 
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(9) 
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— — W„ is ,he minima, possible transit time (hrough gc 
compartment 30 <„ hours,. i, is two-fold lower than the norma, transit time 
Unless the vaiue tha, transit time snoufd approach does ^ ^ 

hour by 1 hour In the direction of the value, [explain??] 

When TPO concentration rises above the threshold, i, sets the value 
which stem cel, transit time approaches, ,o the following number [=*,• 



25 




00) 



where G* is the function 



given in Equation 4, Q is the actual TPO 



concentration and thres is the 



30 



same threshold as for stem cell amplification 

-cover, of committed megaKa oocyte compartments when those [what?, are 

diminished. [?] 
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We assume maxima, possible transit time (300 hours) , above ^ 
transit time does not rise. 

As in amplification rate, in a range of very ,ow oe„ numbers (,00 times 
lower then norma, or ,es S) . the mode, mus, define the sfem oe„ parameter 
^ ,e q uic te s, recover, of SC compartment 30. T h us, transltsc i„ these 
conditions achieves the maximal value, which is 300. 

b. CFU-Meg and MKB compartments: 

Transit times of these CFU-Meg and MKB compartments 40 and 50 are 
completely dependent on TPO concentration. The value, which the transit 
time approaches, is [=?]: 
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where lranslt (cn . UemiKB) m is tne minima| ^ ^ ^ ^ 

Meg compartment 40 (eg. 30 hours, or MKB compartment 50 (66 hours, 
-pectivefy. G is given in Ration 7, and F is the same function as in 
Nation (7 or 8??) , but , he ^ q( ^ ^ ^ ^ ^ ^ ^ 

possible transit time is 100 hours fnr rnn, 

uu nours for CFU-Meg compartment 40 and 350 

hours for MKB compartment 50. 

Ce„ flow between compartments refers no, to the rate of cel, passage 
<rom one compartment to the next, and no, to the number of ce„s tha, 



Cunng a ,ime unit, bu , rather t0 the proportion Qf ^ ^ ^ ^ 

*e next compartment a, any given moment, designated "/W. -Mature" 
cells a,e ones tha, are potential ready t0 pass to the next compartment bu, 
do not necessarily pass. 

> As was noted eartier, i, is assumed that the fraction o, SC compartment 

30 that commits to the megakaryocyte lineage is constant (0.342) and TPO- 
"dependent. From the two Mowing compartments, CPU-Meg compartment 
40 and MKB compartment 50, every mature cell emerges to the next 
compartment regardless o, external circumstances [is mis correct. Or fus, 
■0 TPO circumstances,,. Thus, TPO does no, affec, the « won in the tirs, 
three compartments. 

in con,ras,, ,he fracion of MK16 , MK32, and MK64- megakaryooy.es 
-ha, cominue w„h additional endomitosis and flow to the next compartment 
-ay change from almost o up to t depending on TPO concentration 
Because there is no compartment with pioidy greater than 128N , he 
megakaryocytes ot M K1 28 compartment 90 do no, „ow to any other 
compartment, so there is no devalue for MK128 compartment 90. 

The dependence of MK16, MK32. and MK64 f low0 n parameters on TPO 
concemration is expressed by the following bi-phasic function: 




(12) 



where m.onnorn, is the vaiue of „ OW on under norma, TPO concentration 
(100 pg/ml). This value for MK16 is 0.26, tor MK32 i, is 0.016, and for MK64 



it is 0.4. 



O- * a function ,ha, resembles .he function G in Equation 7. The 

difference is that G° relates fn th= Tm 

relates to the TPO concentration before a certain period 

of «me rather than to the current TPO concentration. The reason for this is 
*a, i, is assumed that a cel, , ha , enters a given megakaryocyte compartment 

5 ' irS ' WhS,her *° "*«- and no, participate 

- Pla,ele, release in this compartment, or whether ,o begin wi, h p, a , e ,e, 
-lease and remain in ,his compartment unti, complete exhaustion TPO- 
dependen, determination o, what fraction of celis wil, ^ each 
occurs a, ,he s,ar, of the cells' path through the given compartment. However 
- the resu,, o, this ^termination can be seen when the ceils tha, "decide- ,o 
'eave ,he compartment aCually leave i,, i.e. after they complete one 
endomi,osis. Thus, ,he fracion o, ce„s ,ha, pass ,o ,he nex, compartment 
«o«on, is calculated based on ,he TPO concentration tha, existed before ,he 
period of one endomitosis (48 hours): 
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(13) 

The time needed for additiona, endomitosis is the same in a„ three 
compartments, MK16, MK32, and MK64. 

III. COMP1 FTP n CT fl„ Fn mr)n 

bone marrow. Each compartment is subdivided into smai, sections that 
the fifth age-section of MKB compartment 50 contains ce„s within MKB 



co~n t 50 ,ha * have co mpartm en, for 5 hour, Every 

hour, a„ the cells in the „one marrow pass t0 the next age . sect| , n , the 

same compartment. 

When ,he ce,l has span, al, the transit time predetermined for i, in a 
» given oompartment, it passes to the hex, oompartmen, to the zero age- 
«*n. Thus, every hour the oeiis tha, leave one oompartmen, fil, the zero 
age-section of ,he next one. The oeiis tha, ieave MK1 2 8 compartment 90 
have nowhere to go and therefore disappear. The zero age-section ot SC 
compartment 30 is .led by a certain traction (0.5,8,) ot the celis tha, ieave 
io SC compartment 30. 

The cells ,ha, release p,a,e,e,s add a certain p,a,ele, number to the zero 
age-section of PL compartment 100 every hour. 

Reference is now made to Fig. 4, which is an illustration 0, the 
.mpiementation 0, the model. The mode, is implemented as a chart of 8 rows 
- and 360 columns. The 8 rows relate ,0 8 eel, compartments, and the columns 
-late ,0 ,he age sections, with the assumption tha, transit time does no, 
exceed 360 hour, This chart is updated hourly according ,0 the rules 
described above. 

representation of the chart of Fig 4 Within th* 

■-■y. 4. within the compartments that 

Proliferation occurs (SC and CPU-Meg,, the number of P ro,ifera,ing cells 
-eases from ,he firs, ,o the las, age-sec,ion. In contrast, the cel, number in 
*e compartments tha, have no proliferating ability remains constant (MKB 
MK128, PL,, or decreases when cells that have undergone additional 
endomitosis leave ,he compartmen, for the next one ,MK,6, MK32 , M K6 4 >. 



Reference is now made to Fig. 6, which is an illustration of another 
representation of the model, based on the time courses of different 
compartments. The rows in the chart represent cell compartments and the 
columns represent time of simulation course. At every time-step of the 
J simulation (one hour of -patient's life"), the number of cells in all age-sections 
is summarized for each compartment and the next column in time-course 
chart (Fig. 6, is filled. Thus, every cell in the chart represents the total 
number of cells in a given compartment at a given time point. 

There is an additional row in the time-course chart that relates to the 
TPO concentration in the blood. TPO concentration is written down every 
time-step concurrently with the cell numbers. 

Reference is now made to Fig. 7, which is a graphical representation of 
the chart of Fig. 6, and is the most useful model output [why? Describe 
graph and how it is used]. 

The implementation of the described model results in a computer 
simulator that describes the changes that occur in the human thrombopoietic 
system (platelet counts, bone marrow precursor numbers, and TPO 

concentration) over several years. The resolution of the simulator output is 

one hour. 

Time units and periods that we will speak about relate to the 
simulated patient's life, not to the running time of the program. 

IV. PARAMETFR.gPl=ri. r| C ADAPTATION 0F MQDFI 

This model may be fit to patients with diverse blood and bone marrow 
parameters. People differ in their baseline platelet counts and numbers of 
bone marrow precursors, in the cells' transit times and amplification rates, 
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rates of platelet release by megakaryocytes, fractions tha, each 
megakaryocyte ploidy Cass donate for additional endocytosis, and in the time 
needed for endomitosi, Furthermore, the rate o, TPO production, receptor- 
and non-receptor-mediated TPO clearance, the threshold of TPO effec, on 
SC compartment 30, and the sensitivity o, different cel, parameters to TPO 
also differ from one person to the next. 

To obtain an idea, fitness of the model to each patient, the paflent- 
related parameters should be g^en individually for each patient. However 
practically, „ would be extreme* difficult to predetermine many of these 
parameters for every patient. Therefore, certain average parameter have 
been calculated based on published data, and are shown in Table , below 
These averaged parameters are used as a framework into which known 
"dividual characteristics are included. Thus, before a particular simulation is 
begun, relevant known information about the individual may be included, 
sometimes replacing certain parameters of the model. 



TABLE 



Feature 



Compartme nt 

SC 



CFU-Meg 



MKB 
MK16 
MK32 



MK64 
MK128 



PL 

TPO-related 
parameters 



Steepness of the 
TPO sensitivity 
curve of 



Baseline 
number 

\3, 



(x 107kg ) 



478 



1250 




5105 
4080_ 
1250 

11 
7.5 



Normal 
amplification rate 
(increase p er honr^ 



Normal transit time 
(hrs) 



1.029 



1.0255 



24 



60 



120 
120 



1 5 050 000 
Rate of 
reduction 

7pg/ml/hr 



SC 

amplifica 
tion rate 



120 



120 



240 

Non - rece PK^ 1 ReTeTtr^ 
mediated clearance Mediated cl^r.n^ 

3%/hr 1 4.776 x 10 10 

f ml/receptor/hour 



CFU-Meg 
amplification 
rate 



CFU-Meg 
transit time 



MKB 
transit time 



Normal cell 
flow* 
(hour -1 ) 



0.342 



1 



0.26 
0.016 



0.4 



Threshold for 
the effect on SC 



2500 pg/ml 



Rate of platelet 
release (hour -1 x 
megakaryocyte" 1 ^ 



_10.2 
20.4 



40.8 



81.6 



Fraction of MK16 Fraction of MK32 
undergoing undergoing 
additional | additional 



Time needed for 

additional 
endomitosis Thrs) 



_48 
^48 
48 



Fraction ofMK64 
undergoing 
additional 



different 
parameters^/) 



30 



100 



endomitosis 



endomitosis T 



«*• - *. log,ri,hm tee is 10, ,ot^ liSfe eC ' ,a °' Sm """ dara ° i » K ampl i ficaUon me <s> is 0.2. 



endomitosis 



Usually, the known patient-related data are no, parameters in ,he form 
def,ned by our mode,, but rather measurements obtained in the clinic (e g 

Cay and vaiue o, post-chemotherapy thrombocytopenia nadir, day and value 
■o o, platele, peak after TpQ adminjstraijon ^ megaka ^ ^ 

Plcdy following some perturbation, etc.,. In these cases, the available data is 
converted into a model-compatible format. 

Sometimes, the only available patient-related data is the graphic 
representation o, the patienfs platele, course following some perturbation 
- (e.g., compressive therapy or TPO administration,. The data may also be 
a P,c.ure of the plateie. course without any externa, disturbance (e.g., cyclic 
•hrombocytopenia,. In these oases the mode, parameters are changed by 
tnal-and-error until a good compliance of the mode, graphic output and the 
Patient, graphs is achieved. „ should be noted, however, tha, even In the 

rather are also based on some analysis [for example?]. 

Specialty, the following tools are available for proving maximum 
flexibility: 

D The user can set the baseline values and a„ other known patient- 

specific thrombotic parameters before starting the simulation. 
2) The user (e.g., physician) can determine how long of a time period to 
simulate, from 12 hours up to several years. 



3) The user can determine .he frequency of showing the course o, a 
patient counts up to the moment. The frequency can change from as 
much as every 12 hours t0 once dun - ng ^ ^ ^ ^ 
simulation. 

> 4) The user can determine the resolution o, the output graph, from the 

houriy representation of the patienfs state down to any other 

resolution. [What does this mean-? r»« .■♦ u 

es mis mean? Can it be resolved at units of 

less than an hour? If so, how?] 

5) The user can choose to view the graphical representation of the age 
distribution through the compartments at any moment of the 
simulation. 

6) The user can imitate a cell-suppressive therapy a, any moment while 
inning the simulation by reducing one or several o, the 
compartments by any value. 

?> The user can simulate exogenous TPO administration a, any 
moment while running the simulation by controlling dose heigh, 
number of dosings or frequency of dosings. 

The simulation tool has been carefully tested with respect to the 

published experimental results anrf h=. 

ai results, and has proved to be well calibrated for 

average human data. P ararneters may be ^ ^ ^ ^ 

efficient use of the system. The following mode, parameters are important for 
individualized adjustment of the model: 

•baseline number of: SC, CFU-Meg, MKB, MK16, MK32, platelets, 
•amplification rate of: SC, CFU-Meg. 
•transit time of: MKB, MK16, MK32, MK64. 



•fraction undergoing additional endomitosis in: MK16, MK32, MK64. 
•rate of platelet release of: MK16, MK32, MK64, MK128. 
•Time needed for additional endomitosis. 
•Rate of endogenous TPO production. 
• Ratio of receptor- and non-recep.or-mediated TPO ciearance. 
•Steepness of the sensitivity curve of: CPU-Meg amplification rate- MKB 
-nsi, time: MK,6, MK32, and M K64 traction undergoing additiona, 
endomitosis. 

[Question: wn.cn o, these parameters can ac,ua„y be measured 
'» *„d input into the program? All? Are there Cher parameters ,00 

such as age, weigh,, medica. history, etc.? , think we „ eed t0 
distinguish between types of parameters.] 

comparison between experimentally obtained data and the simulated mode, 
' Experimentally obtained in v,Vo platele, counts foNowing TPO administration 
are shown in Pig. 8A [ls tnis with ^ ^ 

without TPO is shown in Pig. 9A . Rg , 8B and gB show ^ ^ 

same. By using a TPO schedule designed by the described model, one can 
obtain platele, profiies that are similar ,0 those obtained clinically (Pig 8B , or 
even more effective (Pi, 9B , , n this caS9 , these resu| , s ^ ^ ^ 

administering a pre-calculated TPO protoco, whose total dose amounts ,0 
25% of the original total dose, exptain the po.n, of this better 

The compiete model simuiates cei, and platele, counts in the steady 
state, as we,, as after perturbations ,0 the hematopoietic system, e.g ce„- 
suppressive therapy, recombinant thrombopoietin administration, etc „ is 



possible ,o simulate any protocol of drug administration and any 
hematological state of a patient, regarding his/he, platelet ooun, and number 
of bone marrow monocytes and their precursors. The mode, can be 
adapted to many categories of patients, or healthy platele, donors. „ can also 

one can adjust the mode, so as to yield particular predictions about the 
thrombotic profile o, an individual patient. Platele, disorders, such as 
cyclic thrombocytopenia, may also be simulated. 



EXAMPLE a: Neutrophi, Bone Marrow and Per iph era. Bfood Compartment 
under the Effects o f Growtn-Factors and Trea,men, wft G ranulocyte 
5 Colony Stimulating Factor (G-CSF) 

The neutrophi, ,i„ea 9 e originates in piunpoten, stem ce„s that 
Pirate and become conned ,o the neutrophi, iineage. These oeiis then 
undergo gradua, maturation accompanied w„h probation through the three 
■o n,or P ho,ogica„y distinguishabie mitotic compartments: Mye ,oo,as,s 
promye,ocytes and myelocytes. The myeiocytes then mature and iose the, 
capacity to protiferate, and thus enter the post mitotic compartment. ,„ the 
post-mitotic compartment the ceiis continue their gradua, maturation, which is 
no. accompanied with probation through the three morphoiogicaiiy 
- distinguish^ sub-compartments: Metamyeiccyte, hand and segmented- 
neutrophiis. Ce„s exit the various sub-compartments ,„ the post-mitotic 
compartment and enter the oiood as neutrophil They then migrate trom the 
blood to the tissues. 

The Granuiocyte-Coiony Stimulating Factor (G-CSF) effects an 
> "crease in oioc. neutrophi, levels primari,y by basing product in the 
Mitotic compartment and shortening the transit time of the post-mitotic 
compartment. 

Thus, the „rs, compartment of the mitotic poo, (mye,ob,as.) receives an 
*w o, ce„s from stem-ce,, precursors. ,n„ow for each of the other 

compartments is from outflow of the orevinu* ™* ^ 

8 P reviou s one, subject to multiplication 

factors due to cell replication in the mitotic stages. 



Models regarding granulopoiesis in normal humans and in humans 
with pathologies of the bone marrow were suggested previously in order to 
give a coherent description of the kinetics of granulocytes from experimental 
data (Cartwright - 1964, Mary - 1978, Rubinow - 1974). In recent years 
Schmitz et al. developed a kinetic simulation model for the effects of G-CSF 
on granulopoiesis (Schmitz - 1993), and used i, for the analysis of 
administration of G-CSF to patients suffering from cyclic neutropenia (Schmitz 
- 1995). However, the data Schmitz rests upon for his model has been more 
accurately assessed in recent years by Price et al. (1998) and Chatta et al. 
(1996). Actual empirical data regarding compartment sizes and their transit 
times was not incorporated into their model despite the importance of these 
data (Dancey et al. 1976). 



I- MODEL OF NFDTROPHII , , MC ^ F ftflR 



lEFECTSOFG-CSF 



A. G-CSF 



10 



15 



The effects of G-CSF on the neutrophi, lineage are re|ayed h ^ 

^ " ,hree S * a9eS - The *« * - <*■*■««■ amount of c^ine gi ven 
at time t, which is marked- G' 

adm 

The GMm vector se „ es as , he comrol ^ . ab|e ^ optea(ion ^ ^ 

administration. 

The second stage represents the pharmacokinetic behavior of G-CSF 
« circulation. „ incorporates, tor instance, the haif-iife of G-CSF, and couid in 

activity. This level is marked G' 

blood ' 

G-CSF is eliminated from the blood in a Poissonic manner according to 

WE in "Pharmocokinetics of recombinant human granulocyte-macrophage 
ccony stimulating factor in children after intravenous and subcutaneous 
administration"" Journal of Pharmaceutical Science, 84(7): 824-828 1995- 



G 



Hood =G' blood (\-M2.) + QM 

t y adm 



(14) 



20 



"hare u A is the haHfe of G . cgF ^ b|oo 

blood ^adnt 

Recent data by Terash, K, Oka M, Ohdo S, Furudubo T, ideda C 
Fukuda M, Soda H, Higuchi S and Kohno S, in "Cose association between 
clearance o, recombinant human granulocyte coiony stimulaflng factor (G- 
CSF, and G-CSF receptor on neutrophil in cancer patients', ^microbial 



Aoents and Chemotherapy, 43, 1): 21 . 24 , 1999 , polnts , 0 the - 

of G-CSF effects on the neutrophil lineage, the half-life is considered as a 

5 emerge. 

Only exogenous,,, produced G-CSF is considered t0 affect the Kinetic 
Parameters, and endogenously produced G-CSF ieveis and effects are se, ,o 
™* „ more empirical data regarding the production o, endogenous G-CSF 
- made avaiiable, i, couid be incorporated into the equation as weii. 

The third and flna, stage models the pharmacodynamic effects o, G- 
CSF on the Kinetic parameters. As wil, be eiaborated subsequently the 
dependence or the various Kinetic parameters of the neutrophil lineage on the 

-cave or non-increasing convex function, This reproduces the effects of 
■ oration tha, are seen In clinical studies on the effects o, G-CSF such as 
- study by ouhrsen U, Villeva, JL , Boyd j KannouraKis G, Morstyn G and 
«— D in -Effects of recombinant human granulocyte colony-stimuiating 
factor on hematopoietic progenitor cells in cancer patients", Blood, 7 2(6 , 

2074-2081, 1988. That is addition of C pqp 

s, addition of G-CSF carries a lesser effect when its 

level in circulation is already high. 
B. Biological Model 

Mitotic Cnmp a rf^ nt 

Long-term effects of G-CSF administration taKe place in the 
n«c compartment, though the major contributor to heightened blood 



neutrophil 00unts in the sho „ ^ fe ^ ^ ^ 

shortening of transit time due l0 G-CSF administration, this high ,eve, cannot 

be maintained over the iong term without increased production in the mitotic 



compartment. 



The mitotic compartment is divided into subcompartments. The kth 
subcompartmen, contains aii ceils of chronoiogica, age between W and * 
hours, reiative to the time of entry into the mitotic compartment The number 
of cells in subcompartment k at time f is marked 



as m[ 



(1..T) 



(15) 



where t is the transit time of the entire mitotic compartment, and is assumed 
to be the same and constant for all cells entering the mitotic compartment, 
and /, is a vector reflecting the flow of newly committed cells into the mitotic 
compartment. The biological grounds for this definition is the existence of a 
myeloid stem eel, reservoir, which is known to suppiy new committed ceils to 
.he mitotic compartment. However, the resets actual kinetics are no. very 
well explored empirically. We therefore fix /, to ,eve,s such tha, the overai, 
size of the mitotic compartment as well as the kinetics of the neutrophils in 
circulation would match those obtained empirically. 

Any new bioiogical data that emerges may help define the kinetics 
more accurately within the framework of this mode,, although resuits of this 
model indicate tha, the assumption of a constant rate o, stem cells flowing 



into the mitotic compartment in ,he absence cf G-CSF is piausible. For every 
ne {!.,,, and ,or eve* t amplification occurs a, the exi, from m ' according t0 
Equation 15 as follows: 



<t=m^a n (G' blood ) 

(16) 

where: 

«. is a non-decreasing concave function of G-CSF levels in the blood 
which determines the factor o, amplication in the hourly subcompartmen, n 
«. 'or instance, no ampiification occurs a, subcompartmen, n» a, time X? then 



a =1 



(17)The size o, the morphological sub-compartments in the mitotic 
> compartment at time f is determined as: 

(18)Where n 0 is the firs, hourly sub-compartmen, of a morphological sub- 
compartment and „, i s its las , houriy sub-compartment. The division into the 
morphological sub-compartments is used only for fine-tuning of the Knetic 
parameters with the use of experimental data. 

The mitotic compartment was modeled with an intention to facilitate the 
specific cell-cyce cytotoxic effects of chemotherapy. Therefore, cohorts of 
one hour are modeled as undergoing a process o, maturation and 
amplification culminating in their entry into the post-mitotic as described 
below. Effects of chemotherapy may be incorporated into the model by 
mapping the various cell-cycle phases (G1, S, G2, M) to the hourly cohorts 



mode,ed and ,ormu,a,,n g . funoS , n of the „ ^ ^ 
on these phases. 

The experimental Mature shows wide agreement regarding the 
steady state norma, amounts o, circling neutrohpiis, size o. the p 0 s,-mi,o,io 
> compartment and ,he three morphological, distinct sub-compartments of , he 
n**ic oompartment, and postmitotic transit time and amotion rates in 
.he mitotic sub-compartments (see. tor exampie, Dancey ST, Deubeibeiss KA 
H-tor LA and Finch CA, in "Neutrophi, Kinetics- in Man. Jo uma, o, Clinical 
investigation, 58(3): 705-7,5, 1976; Price JH ^ Gs ^ ^ ^ 

" "E«ec o, recombinant granuiocy.ee coiony-stimuiating facor on neutrophil 
fcnetics in norma, young and eldedy humans', B,ood 88(1): 335-340 ,996- 
and Dresch Mary ,n "Growm traction o, mye,ocytes in norma, human 
granulopoiesis", Ce,, Tissue Kinetics ,9: „- 22 , 1986 , To ^ ^ 
relevant Kinetic parameters, which were either no, avai,ab,e ,n the titerature or 
* were g„en a wide range by experiments, steady state Mnetics was 
-umed and an iterate process was empioyed. These parameters inCude 
.he inflow o, stem ce„s to the mye,ob,as« compartment and the trans,, ,,mes o» 
the mitotic sub-compartments. 

The ha, f „ fe o, b,ood neu,ro P hi,s and the steady state number o, 
neutrophil were taken as 7.6h and 0.4x,0» ce„s* g body weight, respective, 
(Oancey e, a, ,976). S,m,,ar,y, the same ca,cu,a„on may be made tor each 
pa,en, tha, ,s to be mode,ed. This wou,d ai,ow the dynamics o, every patien, 
to be described by the simuiation. The average size o, the post-mitotic 
compartment ( 5,4x,0» ce„s/ kg body weigh, - Da „ C ey, 1976) and the trans,, 
<™ o. ,he compartmen, (1 60 h - Dancey, ,976; Mary, ,966; Price, ,996, are 



compatible with the si>p ann h*if x 

and half-l,fe of the circulating neutrophil compartment 

reported by Dancey, thus supporting the steady state analysis. 

order to determine the amount „, cells in the hourly sub- 
compartments in the mitotic compartment, a„ compartments in the lineage 
5 a steady state assumption. The number ot cells exiting 

'he circulating neutrophil poo, equals the number of cells exiting the post 
n*,ic compartment, which in turn eouals the hourly production ot cells in the 
n-* compartment. Thus, the number ot cells in the ias, hourly cohort ot 
*e mitotic compartment can be de.en.ined trom the neutrophil decay rate 
- which is avai,ab,e in the literature. However, this caicuiation is based on 
assumptions , ha , there is no apoptosis in the post-mitotic compartment 
Direct experimental data by Thiele J, Zirbes TK, Lorenzen , Kvasnic ka HM 

' ^ aP ° P,0SiS ^ PCNA ■ Annals o. Hematology 75(1- 2) . 33 39 

S ~ ^ W °* * - a Phenomenon in normal 

human bone marrow. The size calculated tor the mitotic compartment is 
dose ,0 that experimentally obtained by Dancey and Pnce, thus supporting 

Values tor the production o, ceils in the mitotic compartment can later be 
modified in light of new evidence. 

Regarding the transit time of the mitotic compartment there is iittle 

experimentalists (see Oresch ^ in . Growth ^ q , ^ 
norma, hman granulopoiesis; Cel, Tissue Kinetics t9: n-22, ,986, In order 



10 



15 



•0 determine the tran , t limes of , h6 mWc mon3hoiogfcai 
as m Equation 18, , he following constraints were considered: 

The sizes of the theoretically obtained m orp hologi ca, sub- 
compartments must ,i, , h0S e reported experimentally in norma, 
human hematopoiesis (Dancey, ,976) and under the effects of G- 
CSF (Price, 1996); 

2- A, leas, 24 hours, tne typical cei, cycle, must separate 
amplification points; 

3- The size of ,he las, hourly sub-compartmen, must equal the hourly 
production o, the mitotic compartment (calculated „i, h the 
aforementioned terative process assuming steady state kine,ics,; 

4 ' AmP,i ' i0a,i ° n insWe the °°~n, is se , a, ,he levels 
determined by Mary 1 984; and 

5. The total transit time of the mitotic 
the 90-160 hour range. 



compartment must be within 



20 



above-mentioned constraints. 

» should be noted ,ha, when other alternatives with shorter transit 
«mes were attempted, results couid not be obtained tha. agreed with the 

Furore, a „ be^een our simu(atlon mode ,, s ^ 

"MN counts in P e ripneral blood with emp(rjca| ^ ^ ^ ^ 

achieved without specu.ating extensive* on the nature o, G-CSP effects 



>0 



on no„-co mrnitted stem cel|s (date provided) (( shouw ^ ^ ^ 

"" te ~ Ve «** " available regarding stem oe„, 

The effects o, G-CSF on this compartment are modeled as an Increase 
ra * e °' °* "» — from the uncommitted stem oel, 

5 inCreaSSS ^ - *„d in.roduc.ion of new points o, 

amplification as shown in Equation 15, 16 since little H„ • 

• since little data is available 

regarding the increases in amplification due to G-CSP, an initial a ss ump , ion 
was made tha, amplification reaches ful, potential a, points tha, under norma, 
conditions undergo an amplification of below a factor of 2 . Additional,, i, was 
0 ~ d ^ ^ *« *» - «* sub-compartments and the typical 

the contrary. 

Reference is now made to Figs. X and y, which show a comparison of 

«a.a in the literature. Increased neutrophil production is in accordance with 

the neutrophil counts reported by Price et al m *HH-r u 

Dyrnceetal. In add-on, these increases are 

m accordance with Print's riato ^„ * 

Pnce s data about neutrophil bone marrow pool sizes. 

Reproduction of the effect o, G-CSF on neutrophil counts and the 
n*,ic compartment sizes beyond day 5 of administration was accomplished 
»y assuming an increase „5% with the highest dose o, G-CSF) in the rate of 
cells entering the myeloblast compartment. Alternatively, G-CSF may change 
•he behavior o, the myeloblast compartment such tha, some o, the cells there 
undergo self-renewal instead o, moving on to the promyelocyte compartment 

modified in light o, new experimental data in the future. 



PoSt-mitotin rnm nrtmrnt 



The post mitotic compartment is relatively insensitive to cytotoxic 
chemotherapy. Therefore. I, is biologically accept and computationally 
sensible to mode, this compartment as a single poo, o, cells, such tha, the las, 
hourly cohort o, the mitotic compartment enters the compartment, and a 

proportion of the cells within the compartment enter the neutrophil poo, every 

hour. 

The post mitotic compartment at time t is a single quantity of cells p\ 



10 such that: 



P - h(G' blood ) ■ p' + m ' T 

(19) 

where fc is a convex, non-increasing function of G-CSF levels in the blood 
which takes values in the range of [<MJ. This definition entails 

G-CSF affects the post-mitotic compartment by shortening its transit 
•™ (i. e. decreasing y. Prioe notes that , he ^ Qf ^ ^ ^ ^ 

m «c compartment is no. significantly changed following administration of G- 
CSR This determination is based on counts made on day 5 after G-CSF 
administration. Thus, i, can be safety assumed that any increased predion 
» ot the mitotic compartment flowing into the post mitotic compartment is 
translated over the long-term to an increase in the flow o, cells from the post 
mitotic compartment to the neutrophil pool. This incr9ased „ ow „ 
compensated by increased production in the mitotic compartment only a, a 



later stage. Therefore, an upper limit to the number of 



cells in the post mitotic 



compartment was set, which is at the values 



given as steady state counts (n). 



In brief, the effects of G-CSF on the neutrophil lineage are modeled 
during the first few days primarily as a decrease in the counts of the post- 
mitotic compartment, which is then compensated by an increased production 
in the mitotic pool. This compensation sustains the i 



counts in peripheral blood. 



increase in neutrophil 



10 



5 



Reference is now made to Fig. , 1 , which is a graphical illustration of a 
simulation of the model. Though no empirical date is available on this point, 
simulations of the mode, predict that the number of cells in the post-mitotic 
compartment decreases substantially during the firs, two days of G-CSF 
administration, and then replenishes, so that on the sixth day the counts 
return almost to their normal levels. This replenishment lags behind that of 
Price e, a, report by a few hours. We can thus formulate a testable 
hypothesis, i.e., whether using the same G-CSF protocol Price et al used, 
there is indeed a nadir on day 3 of the treatment. 



Compartment 



Myeloblasts 



Day 0 (no G- 
CSF) (xlO 9 
cells/ kg. 
Body 
weight) 



Promyelocytes 



Myelocytes 



Mitotic Total = 



0.140 



Day 15 of G-CSF (^el^^^ 
treatment (xl(r cells/ 
kg. body weight) 



0.582 



1.373 



Circulating neutrophils 



2.10 



0.4 



0.153 



compartment size 
due to G-CSF 



0.898 



3.564 



4.615 



2.35 



1.09 



1.54 



2.60 



2.20 

5.88 



We shall mark as o'the outflow from the post mitotic 



compartment: 



(20) 

The number o, neutrophils In ,he circulating Wood compartment a, time , is 
marked „' and is modeled as a single quantity of cells, such that: 

n'*' =o' +n'(l-hll) 

(21) 

where h is the half-life of neutrophils in the blood, as defined in the biological 
literature, and is assumed to be constant regardless of G-CSF levels (Lord 
1989), though this could be easily modified. The Kinetics of neutrophils in the 
tissues are not modeled in this work. 



Neutrophils a nn ,„ me Cimila , tnr] p ,„„. 

The elimination of neutrophils from peripheral blood follows a Poisson 
distnbution, and can therefore be described as an exponential function 
(Cartwright, 1 964 , Therefore the rate o, cells leaving this compartment is 
based on half-life determinations available in the literature. Since no direct 
cytotoxic effects of chemotherapy have been described for this compartment 
it is also modeled as a single pool of cells. 

At the normal healthy level we have the following relationship: 

G uood =0,p' =n,m' 

if 



T = m' T =o l = -^xln2 



which reflects the stability of the steady state. 



The Mnetics of G-CSF is also modeled as an exponential distribution 
with a half-life of 3.5 hours (Stute, 1992) (Eq . U ). 

The effects of G-CSF on the Kinetics of the neutrophil lineage appear 
no. to be a linear function of G-CSF administration levels. Since data provided 

5 ' i,era,Ure ^ ,9M > «* ■*» '0 -wo doses (30 and 300 W am/ kg 

body weight, we can only speculate on the effects of other levels o, G-CSF 
After ,ria, and error analysis, i, was found tha, assuming that the effects o, the 
300 W am dose are ,he maxima,, a, the 30-ugram its effects are about 30% 

of the maximum. 

Reference is now made to Figs. 12A and 12 B, which are graphical 
■llustrations o, the effects of G-CSF a. the two doses. The effects as a 
action of G-CSF level are connected piece-wise linear,, This way the 
neutrophil levels observed clinically under both the 300 and the 30-,gram 
protocols are obtained. 



L5 



II. LINEAR IMPI.EMENTATtr.lu n f THE MfinE , 

This mode, will later be Incorporated into an optimization scheme tha, 
*> have as its objective function both the aims of minimizing G-CSF 
adm,nistra.,on and returning the neutrophil lineage to its norma, levels. 

Although the above-ou„ined model may be implemented in any 
number o, optimization methods, linear programming was chosen because of 
* .nneren, advantages compared with other techniques, i.e. its ability to 
prov.de an optima, sCution using partially analytical methods, and therefore 



being more computationally tractable (Gill 199,). On the other hand, 
implementation of mis model in linear programming carries with i, the 
disadvantage that certain computations must be approximated linearly since 
they cannot be performed directly using linear methods. Thus, we shall 
5 compare the 'closeness' of the solution obtained through linear programming 

will be compared with that obtained through another, non-linear method of 

optimization. 

The significant parts of the model that must be modified due to the 
"near programming implementation are the sections in which multiplication of 

o two varies is defined, sinceT^ 

mu.tip,ication is defined as an approximated va.ue constrained within 
Piecewise iinear constraints that most close.y bound the product within a four- 
faced polyhedron in 3-dimensional space whose vertices are 



Where x min , x max , y min , y max are the constant biologica||y defjned mjnjma ^ 
maxima of x and y. 



M (x,y) 



^min X + X rwi ^ ~ 

- ^max X + ^max J ~ 
^min X + x m u y 

- ymax X + *min )> 



X V 

min J min 

X V 
max ./max 



- X 



y 



max J mm 



- X 



mm J max 



(22) 



Multiplication may also be approximated with variations on the linear 
least squares method, by finding one plane thaf s closes, to the four vertices 



defined. 



The other functions that need to be defined linearly are those 
concerning the pharmacodynamics of G-CSF. Due to the nature of these 
functions (either non-increasing convex or non-decreasing concave), these 
effects are implemented as piece-wise linear functions whose breakpoints are 
the doses for which actual experimental data is available (Chatta 1 994). Note 
that the effects of G-CSF on each of the kinetic parameters have not been 
determined in a detailed manner by experimentalists. Rather its effects over 
a few dose levels on the neutrophil blood counts and the size of the 
morphologically different mitotic compartments and the post mitotic 
compartment have been determined. From these data, the effects of G-CSF 
on the actual kinetic parameters (probability of mitosis, transit time and inflow 
of cells into the myeloblast compartment from stem cell progenitors) has been 
reconstructed a. the dose levels available in the literature. These points are 
then connected linearly to obtain piecewise linear functions relating G-CSF 
levels to their effect on those parameters. Further experimental data in the 
future could be used to produce more accurate functions. 

At the amplification points within the mitotic compartment, the linearly 
approximated multiplication operator (Eq 22) is used instead of the product 

defined in Eq. 16. 

At points where no amplification occurs the quantity from one 
compartment is simply transferred to the next according to the following 
Equation: 



*J+1 t 



(23) 

chemotherapy. 

5 THe „ow ou, of tne post mitotic compartmem (£q ^ ^ 

**ned as a linear approximation of a product. 
3- F o rmulation of , he Mode) ^ ^ 0pt|mj2aWon prowem ^ 
Programming 

The simulation spans a finite number of discrete tim* * , 

01 Qlscret e time steps denoted by T. 

10 We d6fine as control variable the vector that 

e vector th at represents G-CSF 

administration a, every ^ ^ 

^ o bj ec,i V e function is defined as nation of , he ,o,lowin 3 expression: 

(Eq. 24) 

5 Where ' " " ™ mber - - * P« «* compartment a, time , and 
' - - -ar we, 9h , in9 CO e„ tefem , The tog ic for ,ormu,a,in g , he o bi ective 

°" ^ iS ^ ^ ~* *° -* - - mitotic compartment, 
-d y state si Z e for a proved penod of time is sufficient for re haM ta,i on 

~ered q uan % o, G .cs, . 

903,8 ^ « — ' — a different we ign , to be 

9 ' V6n Certa ' n ** - » — to Cinicians) , ha , , he , a ,er 
— - ». Post-mitotic compartment sh ou,d be wei gn ,ed more tnan - firs, 



ones. Obviously this is only one of the possible formulations of the objective 
function as defined in the previous section. 

The pharmacokinetics and pharmacodynamics of G-CSF that were 
defined generally in the mathematical model are defined pieoewise linearly. 
5 Some of the considerations that we put into formulating these functions 



were 



based directly on experimental evidence (elaborated in the main body of text). 
We note however, that actual experimental date regarding the direct effects of 
G-CSF on the kinetic parameters in which this model is interested is rather 
scant. Therefore, some formulations were conducted through partly analytic 
and partly trial-and-error methods. 



The formulation of the model in piecewise linear terms will allow use of 
•his mode, as a clinical too, in three ways. Firs,, the mode, will determine the 
effectiveness of various protocols suggested by clinicians prior to their actual 
use on human patients. Second, the mode, allows computation of the optimal 
protocol ,n a given situation of neutrophil counts, so that the neutropenic 
period following chemotherapy is either shortened or completely avoided a, a 
minimal cos. and exposure to G-CSF. Third, the model serves as a 
constituent in a broader framework of clinical tools that win compute the most 
optima, treatment plan for chemotherapy and growth factors. These uses 
should help clinicians administer more rational treatment to their patients 
minimizing both suffering and medical costs. 



Amplification at 
the exit 



2 + 



2 + 



1.5 + 



0 



Mean transit time 
(hours) 



24 



48" 



Size (10 6 cells/kg 
weight) 



0.139* 



48' 



160 



10.96 



0.558* 



1.4' 



5.84 



^l^Knetics under steady state ^a^tt^ humans. 

+ Dresch, 1986. 

-Calculated based on the steady state assumption as elaborated in the main text. 



Compartment 



Myeloblasts 



Promyelocytes 



Myelocytes 



Post mitotic 



Neutrophils 



EXAMPLE 3: Cancer and Treatment with Cytotoxic Drug Delivery 

Introduction 

Cancer is the second leading cause of mortality in the US, resulting in 
approximately 550,000 deaths a year. There has been a significant overall rise 
5 in cancer cases in recent years, attributable to the aging of the population. 
Another contributing factor to the rise in the verifiable number of cases is the 
wider use of screening tests, such as mammography and elevated levels of 
prostate specific antigen (PSA) in the blood. 

Neither better detection nor the natural phenomenon of aging, 
> however, can entirely explain th. increase in new cases of tumors. Meanwhile, 
other cancers, like brain tumors and non-Hodgkin's lymphoma, are becoming 
more common. Their increase could reflect changes in exposures to as yet 
unidentified carcinogens. Current trends suggest that cancer may overtake 
heart disease as the nation's no. 1 killer in the foreseeable future. As gene 
therapy still faces significant hurdles before it becomes an established 
therapeutic strategy, present control of cancer depends entirely on 
chemotherapeutic methods. 

Chemotherapy is treatment with drugs to destroy cancer cells. There 
are more than 50 drugs that are now used to delay or stop the growth of 
cancer. More than a dozen cancers that formerly were fatal are now treatable, 
prolonging patients' lives with chemotherapy. 

Treatment is performed using agents that are widely non-cancer- 
specific, killing cells that have a high proliferation rate. Therefore, in addition to 
the malignant cells, most chemotherapeutic agents also cause severe side- 
effects because of the damage inflicted on normal body cells. Many patients 



develop severe nausea and vomiting, become very tired, and lose their hair 
temporarily. Special drugs are given to alleviate some of these symptoms, 
particularly the nausea and vomiting. Chemotherapeutic drugs are usually 
given in combination with one another or in a particular sequence for a relatively 

5 short time. 

Chemotherapy is a problem involving many interactive nonlinear 
processes which operate on different organizational levels of the biological 
system. It usually involves genomic dynamics, namely, point mutations, gene 
amplification or other changes on the genomic level, which may result in 
increasing virulence of the neoplasia, or in the emergence of drug resistance. 
Chemotherapy may affect many events on the cellular level, such as cell-cyele 
arrest at different checkpoints, cell transition in and out of the proliferation cycle, 
etc. Chemotherapy may also interfere with the function of entire organs, most 
notably, with bone marrow blood production. In recent years molecular biotogy 
and genetics has made an important step (award in documenting many of 
these processes. Yet, for assessing the contribution of specific molecular 
elements to the great variety of disease profiles, experimental biology must be 
provided with tools that allow a formal and systematic analysis of the intricate 
interaction between the genomic, cellular and cell population processes in the 
host and in the disease agent. This system is so complex that there is no 
intuitive way to know how small changes in the drug protocol will affect 
prognosis. But in spite of this intricacy, attempts to improve chemotherapy have 
been carried out by -trial and error" alone, with no formal theory underlying the 
application of specific drug schedules. Such an approach "is apt to result in no 



improvement, only discouragement and liMe US e,u. information for future 
planning" (Skipper, 1986). 



The treatment of oanoer by cytotoxic drug (or drug combination) 
5 deliver is addressed. In this mode,, two generic types of cells are considered: 
host cells and target cells. Targe, cells are, in fact, the tumor. Both types of 
cells may be damaged whi,e exposed to chemotherapy. The aim is to obtain 
the most suitable treatment protocol according to specified conditions and 
Nation, ,, is assumed that the ce,l dynamics are deterministic and known, 
> and that both types of cells are sensflve to chemotherapeutic agents in certain 
Known period of the ceii-cycie phases. I, a cell is exposed to chemotherapy 
dunng par, of * critical phase, there is a chance tha, i, wil, be eliminated 
blocked or affected in any Cher known way. The descnpUon of the dynamics of 
the delivered drugs is assumed to be known as well. 

In order to achieve the goal optimization process is applied to the 
model. The optimization module uses the model predictions in order to search 
for the suitable solution to posted optimization problem. Precise defining o, 
optimization objectives as we,, as the relevant parameters adjustment is done 
according to the settings defined by user/operator for eve„ special case. The 
method can be applied in general cases as well as in specific ones. 



2. Model of Biological System 



The basic layer of the model incorporates a description of age 
distribution of cycling cells and number of resting (quiescent, cells. The term 



■age of the cell- here refers to chio^al age starting from the 
conventional beginning point of mitotic cycle. 

Reference is now made to Fig. 13, which is a schematic illustration of 
the tumor cell cycle layer. The whole cycle is divided into 4 compartments, or 
5 stages (G,, S, G 2 and M). Each compartment is divided into equal 
subcompartments, where |» subcompartment in each stage represents cells 
of age i in the particular stage (i.e. they have spent i time-steps in this stage). 
The quiescent stage is denoted Go. The cell cycle follows a direction as 
shown by arrows (#). Thus, cells enter each stage starting with the first 
) subcompartment, denoted Gi. 

The model can be described mathematically as follows: Let T„ denote 
the maximum duration of t> stage in the cycle. Let . » denote the time 
resolution of the model in discrete time steps. x lm is a function, which 
represents the number of cells in stage k in the I* sub-compartment, at time t 
to U. t. Both time and age are measured in the same unit, in this case, 
hours. Let Q(t) represent the number of resting cells at time , to , + . t. 
Transit) represents the probability that a cell of age / in the stage k will 
move to the next (* + „ compartment. Cells entering the new stage always 
start from the first subcompartment, i.e. from W . This probability may change 
with time, representing the influence of conditions on cycle length distribution. 

By definition, the cell cannot remain more than T k timesteps in the A* 
compartment, as described in the following equation. 

V k : Trans ( k , T „ ) = 1 

The restriction point (R-point) represents a cell's commitment to 
complete the mitotic cycle. Let T R denote the age a, which the cell passes 



through the restriction point in G1. Only cells in d with l<T R can the cycle 
(in the absence of a drug). 

The total number of proliferating cells P(t) can be calculated as follows: 



P(t)= J 

* = G 1 ,S ,G 2 ,M 



1 *! o 



i = 1 



\ 



In every time interval, quiescent cells may return to the proliferation 
pool. Alternatively, proliferating cells may change their state to become 
quiescent if and only if they are in the G1 stage and at age /, where T# i>0. 
To describe this process we introduce the function Guo(i,t) which describes 
the number of G1 cells in age / which become quiescent during time interval 
[t, t+ • t]. This function may receive negative values, accounting for cells that 
return from resting to proliferation. 

As we assume that the exit to quiescence can occur only prior to the R- 
point (even in cancer cells), and that a resting cell that returns to proliferation 
enters the cycle at T 0 , it can be stated: 

V/>r,,W:G^ 0 (/,0 = 0 

It must be noted that this function is not dependent on / and t solely. Its value 
is determined according to current cell distribution and all the general 
parameters that characterize the described cells group. The same should be 
said about the values of Trans vector that can change during the history of 
given population. 



The model traces the development of described group of cancer cells 



using given parameters, by calculating the number of cells in each and 
subcompartment according to the following stepwise equations: 



every 




i 



10 



15 for k=G 2 ,M, k-1 returns the previous stage (e.g. G 2 -7=S). 



20 



25 



These equations make it possible to calculate the number of cells in 
each subcompartment at every time interval A . starting from initial 
distribution (e.g. at time t=0). Since in this model cell ages are measured in 
absolute time units, these measurements refer to the chronological age of the 
cell, and not the biological age, whose units are relative to a maturation rate 
that differs from cell to cell. Consequently, in this model no cell can remain in 
the same age subcompartment after every time step. On the other hand, a 
fraction of the cells that leaves any subcompartment may be transferred to 
the first subcompartment of the next stage, according to probability vector 



Trans(k,i,t). This vector provides thp ah.iih/ ^ 

yioviaes me ability to account for variability of cycle 

lengths while retaining a deterministic approach. 

The behavior of the cell population in this mode, is completely 
controlled by two components: Trans vector, and 01 *G0 function. These two 
5 functions determine uniquely the outcome of every single ,i me step and 
conseouen.lv the result over long periods. Thus, ,hey are referred to as 
■contro, funcons-. The values o, these functions may be dependent no, only 
on time and age of cells, but also on the current population state (or 
generally, on the whole history of the population, as we,, as on the 
environment associated with a given cel, group. However, those parameters 

here is suitable for describing highly homogenous group of ce„s. Therefore 
.he basic ,ayer of the mode, should give a realistic description for a uniform 
group of cancer ce„s for which environments, conditions and relevant 
biological properties are defined, in a way that wil, a„ow the construction o, 
the control functions for the group. 



• The General Tumor Model 



in the genera, approach the who,e model is viewed as constructed from 
amilar components, each of them derived from the basic structure described 
in the previous paragraph. Each component represents ce„s that are 
subjected to the same environmental conditions and, therefore behaves 
similarly (to be denoted homogenous group). The whole tumor is modeled 
as a union of many varying homogenous groups of cel,s, where the 



10 



development of each 9r0U p can be accurate* priced (when loca, 
conditions are known). 

This genera, mode, simu,a,es progress of the tumor in discrete steps of 
*me. At each step the number of ce,,s in each subcompartmen, of each 
group is ca,cu,a,ed according to the previous state, parameters of tumor drug 
concentration, etc. The parameters of the tumor must inciude a„ the 
information that is reievan, to prognosis. Some o, these parameters are 
defined ,oca„ y , e.g., those reiating to the tumor's geometry. For this reason 
the representation of the spatial structure will be included. 

The cells will be able to pass between the groups during the 
development o, the tumo, This a„ows the representation of the changes in 
*e loca, conditions during the tumor evolution (e.g. forming o, necrotic core 
improvement in "living conditions- in vascuiarized regions, etc). In addition 
all the parameters of the tumor may change in accordance with the dynamics 
of the cancer. 

The calculation of the tumor development over time wil, be done by 
stepwise execution o, the described simu.ation and can be used to predict the 
outcome of the treatment or in fitness function for search algorithms 



2.3. From General to Individual Tumor Model 



When the genera, theoretical description of the mode, is accompiished 
by adjusting a„ the parameters that determine the behavior o, the modeled 



th. task we will estabiish the connections between mathematical parameters 
(most o. them will have direct biological implication) and every kind of data 
.ha, is practically obtainable in the clinic. These connections may be defined 
5 through research on statistic correlation between different parameters 
("eluding genotype-phenotype correlation,, or using advanced biochemical 
research (which may establish quantitative relations). 

Thus defined, fhe mode, wil, be able to give realistic predictions for 
treatment outcomes either for specific patient or for a broad range profiles o, 
» patients and diseases. This too, can serve to perform the prognosis of either 
an untreated cancer patient, or as a basis for treatment modeling as is 
described below. 



3. Introducing Pharmacology 



in order to simulate cancer treatment we add to the above model the 
pharmacologic component. We mode, pharmacokinetics as well as 
Pharmacodynamics for specific anticancer drug, We begin by representing 
cell-cycie specific drugs, however the model implies no restriction on the type 
of drugs to be employed. 

We model the distribution o, the drug in and around the tumor as well 
as in the blood (the drug kinetics,. For this purpose, we use the suitable 
mode,, defining i, precisely for every certain type of the drug. The 
concentrations of drug in the body are calculated a, every time step in 



accordance with the drug administration specified by the protocol, and 
different processes that define drug kinetics in the body. 

The dynamics of the drug are represented through the direct influence 
of the drug on tumor cells. The effects on the proliferating cells are mostly 
5 blocking the cycle in different stages (which can be modeled as cell airest) 
and cell death (immediate or after being in the block). Cell-cycle specific 
drugs are believed to have no direct influence on quiescent cells, but can 
affect them indirectly by killing proliferative cells and therefore changing local 
conditions. Where additional types of drugs added to the model, their effect 
> on any kind of cells may be ,00 modeled as killing certain fraction of cells 
(which is dose-dependent) or changing the behavior of the cells. 

Additional phenomena that may prove significant in drug kinetics and 
dynamics (e.g. rate of absorption by the cells, development of tumor 
resistance to specific drug, etc.) can be introduced into the model ,0 make it 
as realistic as needed. 

The description of the drug in the model is done in terms of quantitative 
functions, which enable to calculate the drug amounts at certain locations and 
the tumor response ,0 i, at every time step. In the general case, these 
functions include parameters that depend on the specific data (drug type, 
body parameters, characteristic of the tumor, etc.) and can be determined in 
given situation (patient/case). The relation between clinical data and these 
parameters can be established in the ways similar to those described for the 
cancer model. 

The combination of cancer model with the drug model described above 
makes i, possible to predict the outcome of the treatment, given the relevant 



parameters ,or the drug, ,he cancer and ,he patient Again, the prognosis may 
be made for specified cases as weli as for broad profiles of patients or 
disease. This simulation also serves to build the fitness function used for the 
optimaization objectives. 



4. Combining with minimising hnst tnyi-ity 



Although an accurate predictive tool, the model that represents 
> chemotherapy of tumor alone cannot be used in optimization, for i, posts no 
constrain, on the choosing o, the treatment. Actually, this model impiies using 
as much drug as possible until the final elimination of the tumor; while in the 
live system the toxicity of the drug is the most important constrain, limiting the 
treatment. In mos, cases o, anticancer chemotherapy the dose-limiting 
toxicity is bone marrow suppression, the two mos, sensitive bone marrow 
lineages being granulopoiesis and thrombopoiesis. Accordingly, ,hose two 
were chosen as an example and are modeled separately and in a similar way 
'o predict the negative effect of the chemotherapy on them. These models 
reconstruct the damage caused by the chemotherapy to the bone marrow 
cells and the recovery of these lineages (treated by specific growlh factors,. 

Thus, the whole system is capable of predicting the result of 
chemotherapy treatment for the tumor as well as for bone marrow cells, 
allowing the use of the protocols that combine anticancer drugs and growth 
factors for healthy cells. 



Chemotherapy toxicity ,o any other normal host eel, population can be 
similarly taken into account, i, „ is defined as relevant tor dose and schedule 
optimization. 



VI. INDIVIDIIAI I7ATIQM QF THF Mnn F | ^ 

Due to a great degree of heterogeneity between malignant tumors (even 
among simiiar tissue types) and between patients, it wou.d be advantageous 
to adjust the treatment protoco. to the individual case. This individualization 
10 procedure includes three aspects: 

1) individual parameters of tumor dynamics 

2) individual parameters of patient-specific drug pharmacokinetics 

3) individual parameters of the dynamics of dose-limiting norma, host 
tissues. 

Relevant data concerning individual cases can be obtained by research 
on statistic correlation between different parameters (including genotype- 
Phenotype correlation), or using advanced biochemical research (which may 
establish quantitative relations,. In the genera, mode,, important dynamfc 
parameters are estimated from experiment studies conducted in certain 
» patient populations. Any of these parameters, when available on the per 
patient basis, can be individuated, white those that are unavaiiabte can be 
'eft as a P opu,a.ion-based figure. This approach allows continuous increase 
i" the degree of individualization of the treatment protocol with progress in 
the technology of parameter evaluation (e.g., oncochips). 



All different parameters may then be adjusted, which will result in an 
adjusted array of models to be simulated. Parameters may include many 
different factors, which are adjustable according to the needs of a 
Pharmaceutical company for general use of the treatment, or may even be 
5 individualized for use by a specific clinician for a particular individual. Examples 
of parameters may include, bu, are no. limited to age, weight, gender, previous 
reaction to treatment, desired percentage of healthy body cells, desired length 
of treatment protocol, pathologic or cytologic specifics, molecular markers, 
genetic markers etc. 

In order for the system to be user-friendly, all possible parameters are 
termed in ways that are easy for the user to understand. 

Once all the parameters are set, an array of solutions is produced 
based on the input parameters. A number of possible protocols can be set (is 
thus generated by the computer), a fitness tacton is applied, w^ich results in 
scores for each of the proposed solutions. 




Refer back to Fig. 2. This model makes it possible to check any given 
treatment protocol and to choose a very goc<, one according to use.s criteria. 
The user may be a physician, a drug developer, a scientist, or anyone else 
who may need to determine a treatment protocol for a drug. The specific 
parameters may include several categories, such as individual patient 
characteristics and/or medical history, needs of a specific user (research. 



efficacy, treatment, etc.), and other particulars (such as maximum length of 



treatment, confidence level, etc.). 



■ That is, an array of possible treatment protocols is created from which 
the optimal treatment protocol can then be chosen. I, should be noted that the 
method does not imply the fitness estimation for all possible protocols. The 
use of operation research allows a much more sophisticated, ye, resource 
saving procedure. 



cancer 



An example of this procedure will be described as it relates to 
treatment by chemotherapy, as described in the third embodiment of the 
invention above. However, i, should be noted that a similar procedure may be 
performed in any of the embodiments. 

The procedure implements cell growth and cell death procedures, as 
defined in the detailed model above. There are certain pre-defined 
parameters, including the lengths of the host and target life-cydes, the 
lengths of their critical phases, and a resolution factor, tha, determines the 
length of a single time unit. The user is asked to define an action (treatment 
or non-treatment). Simulation of cell growth and death is then performed for 
the single time unit This procedure is repeated until the end of the total 
simulation. Alternative*, the choice of treatment or non-treatmen, is made by 
the processor, with many possible permutations considered. In that case, the 
protocol space would be ve,y large, and the resolution would depend on the 
(selected) length of the time unit (and computer capacity) 



There are two procedures: one for growth simulation during treatment 
and one without. The array in which the numbers of cells are kept is updated 
once per time unit, whether with or without treatment present at that time. 



VII. Defining the Fitn«. s Function 

The fitness function is an important tool in Operation Research. In this 
case of protocol optimization, it allows the comparison between a number of 
different protocols each one of them scoring differently with respect to various 
objectives that can be set by the developers or by the users and identifying 
the protocol for which the highest weighted score is predicted. The fitness 
function calculates for any given protocol its relative efficacy ("score" or 
"fitness"), thus enabling a definite decision of the best protocol from a given 
set of protocols. 

In different cases, different objectives can be formulated. There are 
several settings in which such a model can be used, including bu, not limited 

to: 

One) clinical practice- where objectives can change depending on 

type of disease, condition of the patient, purpose of treatment, etc. 
Two) pharmaceutical company- where objectives can be aimed at 

finding the therapeutic window and an optimal schedule. 
Three)scientific setting- research oriented objectives can be aimed at. 
In each case, a particular fitness function can be formulated, reflecting 
all given requirements. Thus, in any particular case one can compare 



between different protocols and obtain the most suitable to his/her special 
purposes and needs. 

Examples for some alternative objectives are given: 

1 ■ The smallest number of drug dosings required for achieving any 

5 given aim. 

2. The lowest total drug dose required for achieving any given aim. 

3. The minimal total amount of drug needed for rehabilitation. 

4. The smallest deviation from the baseline at normal cell 
population count (e.g., platlet nadir) after chemotherapy or another cell- 

> suppressive treatment. 

5. The shortest period of disease (e. g-.thrombocytopenia). 
Using the fitness function it is possible to a) estimate the efficacy of a 

given protocol, b) search for the solution of an optimization problem, i.e., 
predict which protocol will be best of many potential protocols considered for 
curing/relieving the patient. 



VIII. Solving the Qptimi^tion Problem 

The optimization problem is stated using the described models: to find 
the protocol for drug administration (with option to growth factor 
administration) which maximizes the given fitness function. 

As explained above, we define the fitness function according to the 
user requirements. For example, the goal of the treatment may be defined as 
minimizing the number of cancer ce.ls at the end of the treatment, minimizing 



the damage to the BM cel.s throughout the treatment or at its end, and curing 
the patient (where cure is defined precisely) as quickly as possible. Note that 
the fitness function may also include goals such as maximizing life 
expectancy, minimizing cost of treatment, minimizing treatment hazards 
5 and/or discomfort etc. Generally, the aim of optimization is to find the best 
protocol, i.e. the protocol that generates the best value of fitness. 
Customarily, this is achieved by mathematical analysis. However, 
mathematical analysis is restricted to over simplified models, whereas, in 
order to accommodate biologically realistic parameters, the described models 
are very complex and, therefore, cannot be solved analytically. On the other 
hand, the practical purpose of the treatment is not to find the best possible 
protocol (i.e., the global optimum) but "only one that will suit the user's 
objectives, even if its fitness is not absolutely the best (i.e., the local 
optimum). For this reason we can be satisfied with the solution that can be 
shown to promise the pre-defined objectives. 

Hence, the optimization problem may be reformulated as follows: for 
given initial conditions, find the treatment protocol which will fulfill the user's 
requirements (e.g. curing a patient according to given definitions of cure) and 
subjected to given limitations (e.g. treatment duration, drug amounts, etc.). To 
this end it is not compulsory to find the global solution. It is enough, with 
regards to objectives and limitations, to perform search, using search 
algorithms, in certain regions of the protocols' space, and find the local 
maxima of the fitness function. After determining the locally best protocols, we 
can verify that they serve one's objectives and check them numerically for 
stability. 



Such a strategy can be used for identifying patient-specific treatment, 
as well as in the general case, where only the profile of the disease and the 
drug are specified. If more patient-specific data are supplied, the solution will 
be tailored more specifically. On the other hand, the optimization program 
could propose general recommendations for the protocol types for certain 
kinds of disease, treated with a certain kind of medication. 

It will be appreciated that the present invention is not limited by what 
has been described hereinabove and that numerous modifications, all of which 
fall within the scope of the present invention, exist. For example, while the 
present invention has been described with respect to certain specific cell 
lineages, the concept can be extended to any other lineage and treatment 
protocol which can be detailed mathematically (e.g., viral or bacterial diseases). 
Furthermore, certain assumptions were necessarily used in computing the 
mathematical models of the embodiments. Values and equations based on 
these assumptions can be changed if new information becomes available. 

It will be appreciated by persons skilled in the art that the present 
invention is not limited by what has been particularly shown and described 

herein above. Rather the scope of the invention is defined by the claims which 
follow: 



What is claimed is: [independent claims listed now, dependent claims to 
follow ] (Note that we are not in a p osition to revise the claims. 2. Ag iir) 

1 . A system for individualized optimization of a treatment protocol, the system 
comprising: 

a system model comprising: 

a model of a biological process; and 

a mathematical model of treatment effects on said biological 

process; 

a plurality of treatment protocols; 

means for adding individualized parameters to said system model, 

whereby said system model is modified based on said individualized 
parameters; and 

a selector for selecting an optimal treatment protocol from said 
plurality of treatment protocols based on modification of said system model. 

2. A system for identification of optimal treatment protocols using heuristic 
analysis, the system comprising: 

a biological model; 

a plurality of treatment protocols included within said model, thereby 
providing a plurality of treatment models, wherein said treatment models 
provide effects of said treatment protocols on said biological model; and 

heuristic means for identification of said optimal treatment protocol 
from said treatment models. 



Wha, is Caimed is: tindepend9nt c)ajms (jsted _ depen(jent 

' 0W1 <M *^^ 

1 • A system for individualized optimization „f a ♦ . 

comprising: ** protoco '. *» system 

a system model comprising: 

a model ofa biological process; and 



process; 



a mathematics mode, of treatment effects on said biological 
a Plurality of treatment protocols; 
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-ans for adding indlVWuato parametere ^ ^ ^ ^ 

whereby said 5ys ,em mode, is modified ^ 
parameters; and 



on said individualized 
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' rr: rr: 

a biological model; 



3. A system for presenting a individualized 
system comprising: 



model of a biological process, the 



a processor for generating a model of a biological process; 
a Plurality of experimental date, included within said model; and 

means for including individual data into said processor, thereby 
Proving sa W indicated mode, o, sa W biological process. 



4. A system for prediction of outcomes of a treatment 

comprising: 



protocol, the system 



o 5. A system for describing treatment effects on a bfologica, process 

6. A method for individualized opflmizafen of treatment, the method 
comprising the steps of: 

providing a biological model; 
inputting parameters into said model; 

computing treatment protocols for said mode, within said parameters; 
applying a fitness function to said treatment protocols; and 



choosing an optimal treatment protocol based 



on said fitness function, 



7 thetptoT ° PWZa,ton * ^ - ™« ~, 

administereT" 9 ' " " ^ « 



providing behavior characteristics of said administered 



medication; 



of a ~ iine ~ 9 — — ■ 



calculating a 
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determining a treatment protocol 
M l,n 1 1 r °J> th e system comprising: 

a system model comprising: 



a mathematical model of a neutrophil lineage; 



and 



a mathematical model of efferfe nf 
lineage; ° f treatmen * on said neutrophil 

a plurality of treatment protocols; 

means for adding ind , iduafeed ^ ^ ^ ^ 

P ^ ° f ~ S » ~on of said system model 

11. 
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Simulations showing that if the protocol is pre-calculated then a similar or a 
higher efficacy can be obtained using 4-fold reduced total dose of TPO. 
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Fig. 8: TPO given to healthy donors- Results of TPO clinical trials from recent research on healthy platelet donors as 
25 compared to our computer simulation results. Arrows indicate the start of TPO treatment. (A) Comparison of experimental 
data from published articles 1 (black) and our model simulation (green), in both TPO was given as a single IV dose of 1 2 
ug/kg on day 0. (B) Comparison of the same experimental data (black) and our proposed TPO administration protocol- the 
total dose in the simulated protocol was 0.3 (ig/kg (blue). 
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Abstract [return to Contents ] 

Thrombocytopenia is a common hazardous blood condition, which may 
appear in different clinical situations, including cancer chemotherapy, where it 
represents a major dose limiting side effect For this reason, and due to the 
constant demand of platelet harvesting for transfusions, a specific control of 
blood platelet counts is strongly needed. 

Recently, such a thrombopoiesis-controlling cytokine, thrombopoietin 
(TPO), was isolated and its human recombinant analog is now available. In 
order to accelerate a search for the most effective TPO treatment protocol, and 
to reduce the investigation costs (including number of potentially endangered 
patients), we offer a mathematical model, which is implemented as a computer 
simulation. Our model reproduces the dynamics of the human thrombopoietic 
lineage in the bone marrow, of platelet counts in the periphery, and effects of 
TPO administration on them. 

Using this model, we are able to develop general and patient-specific 
optimal protocols of TPO administration for minimizing post-chemotherapy 
thrombocytopenia and for maximizing platelet harvest in platelet donations, 
while reducing TPO treatment cost. 

Introduction [return to Contents ] 

Aim [return to Contents] 

The most harmful dose-limiting side effect of cancer chemotherapy is the severe 
suppression of the hemopoietic system. This results in a state of immune deficiency 
and reduced blood coagulation as granulocyte and platelet counts fall drastically. 1, 2 
Therefore, post-chemotherapy rehabilitation of hemopoietic system is of great 
importance. However, most of the growth factors that are currently used (e.g., EPO, 
G-CSF) are almost impotent in stimulation of thrombopoiesis. This led to an intense 
search for specific thrombopoiesis-stimulating factors in the previous decades. This 
search was motivated also by the constant demand for platelet donations for 
transfusions to thrombocytopenic patients. 

Recently, such thrombopoiesis stimulating cytokine, thrombopoietin (TPO), 
was isolated, and its synthetic analogs became available, and are under clinical 
investigation. 2, 3f 4 ' 5 ' 6 These compounds were shown to have the same biological 
activity as TPO, thus the term hence the 

The term TPO will be used without distinguishing between its different forms 
and analogs. It was shown to have a wide spectrum of useful effects, almost without 
adverse ones. 

Most authors that studied TPO effects on thrombopoiesis emphasize the critical 
role of the TPO administration protocol in chemotherapy-associated treatment and in 
platelet harvesting. 2, 4 ' 7 ' 8 ' % 10 However, the data accessible for selecting the proper 
protocol is very insufficient, because of TPO's recent discovery. Therefore, it is of 
special interest and importance to model thrombopoiesis and TPO effects on it 
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mathematically. An implementation of such a model as software may accelerate a 
search for the most effective TPO treatment protocol. 

In this work we build such a mathematical model and implement it as a 
computer simulation that imitates a course of patient's bone marrow progression 
(particularly, compartments relevant for platelet production), peripheral platelet 
counts, and TPO concentration changes in the system for a period of several years. 
Both cell-suppressive treatment (chemotherapy, for example) and exogenous TPO 
administration, are considered in our model. The model allows their introduction 
concurrently, in any possible combination and timing. 

Currently, our model enables a search for an optimal protocol by trial-and-error. 
An additional module is to be implemented shortly, which will use Operation 
Research tools while running the simulations, to find the best protocol. The great 
contribution to the traditional biological and clinical experiments is its cost and time 
reduction, and the decrease in need for experiments on animals and/or humans. 
Additionally, the optimal protocol search may be accelerated by some analysis of 
model tendencies. 

Modeled biological System [return to Contents] 
Thrombopoiesis [return to Contents] 

[ REFERENCE to a textbook of homo/thrombopoiesis (to look for thorn in (50)) 
Review of thrombopoiesis see in (35Y additional data see in (11. 12. 13, 3L 31-Kuter 
et al (1997V): brief review see below G 1): Branohog ot nl (1 075)(?) J Kntnr nt nl n 007Y 

] 

Like all other hemopoietic lines, the thrombopoietic line originates from poorly 
differentiated, multipotential cells, that are capable of some division and self- 

* * 3 1 3 133 

reconstitution. — 1 Such bone marrow cell compartments as pluripotential | 

hemopoietic stem cells and common myeloid progenitor cells (CFU-GEMM) have 
more or less these characteristics.— 

Gradually, the cells become more and more differentiated and thus committed 
to the thrombopoietic line. At this stage they proliferate extensively. Colony-forming 
Units - megakaryocytes (CFU-Meg) is an example of such compartment. Sometimes 
fee-burst-forming units - megakaryocyte (BFU-Meg) . nromegakarvoblasts or the 
megakaryoblasts^ compartments are considered as having the same similar 
properties. 11 ^ 

The committed thrombopoiotic m egakarvoevtopoietic cells, megakaryocyte 
precursors, go through several stages of maturation when the proliferation is 
impossible. However, megakaryocyte maturation is somewhat different from that of 
other hemopoietic lines. Here, along with cytoplasmic maturation, cell nuclei undergo 
mitotic events. However, although the DNA material of these cells doubles, cell 
division does not happen. Such incomplete mitosis is termed endomitosis or 
endoreduplication . Consequently, the cell becomes poliploid with 2N, 4N, 8N, etc., 
amount of DNA. Some authors call the cells with 2N to 4N chromosome number 
promegakarvoblasts. others call them megakaryoblasts orr oth e rs call thorn immature 
megakaryocytes* 4 ". 11 ^ 
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Usually, megakaryocytes do not start to release platelets until they reach the8N 
to 16N state. 1 Then they begin to create demarcation membranes that envelop 
cytoplasm fragments generating platelets. The platelets are released into the blood 
stream. A small fraction of the megakaryocytes do not cease their endoreduplication 
at the 16N-stage, but rather continue with one or more additional endomitoses and set 
thus a ploidy of 32N or more. — The amounts of cytoplasm, cell volume and the 
ability to release platelets increase proportionally to the cell ploidy. 12, 13 

Thrombopoietin [return to Contents] 

The major sites of TPO production are the liver and the kidneys. TPO is also 
produced in the spleen and the bone-marrow, but the production rate in these organs is 
significantly lower (5 times in mice) than in the former (???). 2, 14, 15 ' 16 ^ rSome 
authors also find low TPO production in many other sites in the body. 15xi2 The rate of 
liver and kidney TPO production is constant under thrombocytopenia and 
thrombocytosis of every examined severity. 2, 14 ' 15 ' 16 Spleen and bone-marrow 
production, on the other hand, is inversely related to the megakaryocyte and platelet 
mass. 2, 14 However, these changes are extremely small compared with total TPO 
production, and are therefore considered negligible. The possibility of high local TPO 
concentration in these organs, in which megakaryopoiesis occurs, was rejected by 
experiments that showed that TPO concentration in the bone marrow equals that in 
the peripheral blood. 17 Normal TPO concentration varies significantly between 
healthy subjects. 34 However, intraindividual variation in steady state TPO 
concentration is relatively low. 34 

TPO-receptors on the platelet and megakaryocyte surfaces are the main TPO- 
clearance mechanism, and thus their major regulation mechanism. Hence, TPO 
concentration is inversely related to the total platelet and megakaryocyte mass. 2, 5 ' 17, 

The effects of TPO on the thrombopoietic line can be divided into three kinds: 
(i) stimulation of proliferation of the megakaryocyte progenitors which have an ability 
to proliferate; 7, 9 * 2l > 22 > 23 - 3I - 32 - 33 (ii) stimulation of maturation of megakaryocyte 
progenitors; 7, 22 , 31, 321 33 (iii) induction of an additional endomitosis of already mature 
megakaryocytes that leads to the increase in the average megakaryocyte ploidy 3, 21> 22, 
23, 31, 33 . These effects are mediated through the same TPO receptor. 19 However, the 
mechanisms underlying these effects may be different. Hence, it seems possible that 
one of them may be introduced when another is saturated, at least partially. 

Methods [return to Contents ] 
Model Structure [return to Contents ] 

The search for the best protocol of TPO administration was divided into two 
major parts: 1) creation of the model and computer simulation, which closely imitate 
thrombopoiesis of any given individual under different circumstances detailed below; 
2) producing the module that can generate the best protocol for general purposes, as- 
well as for the specific patient, according to a criterion specified by a physician, a 
scientist, a pharmacologist, or any other user. The optimization scheme uses a flexible 
fitness function, and is based on the results of simulation supplied by the first part. 
This paper will discuss only the first part. 
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First, we had to generate a scheme of thrombopoiesis that can be expressed 
mathematically. This scheme has to include certain biological mechanisms and 
parameters in order to be able to give a realistic prediction (Fig. 1 - ???). 

Scheme of the Thrombopoiesis Model (return to Contents ] 

We divided the entire thrombopoietic lineage into eight compartments. The 
first, called abstractly "Stem Cells" (SC), refers to all bone marrow hemopoietic 
progenitors that have an ability to differentiate into more than one line (e.g., 
pluripotential stem cells, CFU-GEMM, etc). Cells of this compartment in our model 
proliferate, mature, and then differentiate into megakaryocyte progenitors or other 
precursors, or give rise to "new" stem cells. Although the consideration of the last 
process, i.e. the renewal of the stem cells by "new" ones, is not completely correct 
biologically^^ 2 , it may serve as an_acceptable assumptio n since we do not elaborate 
the characteristics of this population in details . 

Cell death through apoptosis may have a significant effect on cell numbers in 
the proliferating compartments.— However, we do not consider this process 
separately, but rather include its effect together with the effect of cell proliferation 
into total amplification of cell number ( s c) in a specific compartment. As far as we 
know, there is no data about significant apoptosis in the non-proliferating 
megakaryocyte compartments. We therefore disregard apoptosis in these 
compartments. 

Kinetic data regarding the various sub-populations of the stem cell 
compartment (i.e. rates of proliferation, maturation and self-renewal) is rather scarce. 
We therefore assume one homogeneous stem cell compartment with biologically 
acceptable parameters. Based on some works that showed that probabilities of stem 
cell differentiation into one or another hemopoietic lineage are constant, 24, 25 ' 26i ^ we 
assumed that a proportion of mature stem cells that flow into the megakaryocyte 
lineage ( sc) is fixed. The same was assumed about stem cell self-renewal. Thus, 
after the cells spend a defined transit time in the Stem Cell compartment (sc), 
certain constant fraction of them returns to their "young state", i.e. starts their passage 
through the SC compartment again. Another fraction ( sc , also constant) passes into 
the next compartment ealted-nani^"CFU-Meg". The remaining "stem cells" are 
supposed to go to hematopoietic lineages other than megakaryocyte. 

By the name CFU-Meg (Colony-Forming Units - Megakaryocytes), we refer to 
all the cells that are already committed to the megakaryocyte line but are still capable 
of proliferation. Cells of this compartment in our model, like those of the previous 
one, spend some time multiplying at their own amplification rate cf u and maturing 
before they lose their proliferative ability and pass to the megakaryoblasts (MKB) 
compartment. 27 The time they spent in CFU-Meg compartment is C f u . 

The MKB compartment in our model includes all the cells that have lost their 
ability to proliferate, but are not sufficiently mature to release platelets. We assume 
that megakaryocytes do not start to release platelets until they reach the 16N-ploidy 
phase. 11 Thus, "megakaryoblasts" is more of a conceptual name; it refers to the 2N, 
4N and 8N cells of the megakaryocyte lineage that cannot divide, at all stages of 
cytoplasmic maturity. 16N-megakaryocytes, before they start to release platelets, are 
also included in this compartment. After these cells spend their designated transit time 



( mkb) in the MKB compartment, they move to the next bone marrow compartment, 
termed MK1 6. 

The cells of the MK16 compartment are megakaryocytes of 16N-ploidy class. 
They release platelets at a constant uniform rate (r M Ki6) until they exhaust their 
capacity (c M ki6) and then are disintegrated. 12 '~^ 27 We assume that the € cell volume is 
in linear relation with megakaryocyte ploid v based on some evidences from 
literature. 12 Hence, we assume that all 16N-megakaryocytes have the same volume at 
the start o f their path through the compartment and, thus, the same platelet-releasing 
capacity (cmkks)- Therefore, we suggest that all platelet-releasing 16N- 
megakaryocytes spend the same transit time ( M k) until they are exhausted and 
disintegrated. 

However, not all 16N-megakaryocytes are destined to this fate in our model. 
Part of them ( mkis) does not participate in platelet release, but continue with another 
endomitosis and became 32N-megakaryocytes. These constitute a distinct 
compartment, MK32. Hence, after the cells spend the period needed for additional 
endomitosis ( ) in MK16 compartment, a certain fraction of them leave it and pass to 
the next one. 

32N-megakaryocytes also release platelets. The rate of platelet release is 
constant for every compartment and proportional to the ploidy-state of the 
megakaryocytes in it. 12, 13 Thus, 32N-megakaryocytes release platelets at the rate 
rMK32, which is twice that of the 16N-megakaryocyte. However, we suggest that 32N- 
megakaryocytes are not exhausted sooner than the 16N class since they have twice the 
volume and platelet-releasing capacity (c M K32). Consequently, we assume that all 
platelet-releasing megakaryocyte compartments have the same transit time (mk)* 

Again, a fraction of the cells ( M K$i) is not engaged in platelet formation, but 
continues to the 64N-stage. We assume that the additional endomitosis takes the same 
time ( ) as in the MK16 compartment. 64N-megakaryocytes continue this behavior, 
and a certain fraction of them ( M K64) becomes 128N-cells.— ^ The time needed for 
an additional endomitosis remains the same () . Normally, 6 4 N — and 128N 
megakaryocytes ore not found in humans. [REFERENCE] Howovor, thooo colls were 
d e tected in human bono marrow in some abnormal conditions. [REFERENCE] 
Therefore, wo includ e d those compartments in our model, but in the normal otato the 
number of colls — in thorn is negligible. — According to current knowledge, 
megakaryocytes of greater ploidy classes were not encountered in human. 

Finally, a platelet (PL) compartment is introduced into our model. This is not a 
bone marrow compartment, but rather the platelet pool in the peripheral blood. 
Platelets, which are released from megakaryocytes of 16N-, 32N-, 64N-, and 128N- 
ploidy classes, enter this compartment. There are two mechanisms of platelet 
elimination from circulation: By age-dependent destruction and by the normal 
utilization diroctod at (???) in order to maintaining the integrity of blood vessels. 23, 28 | 
The first mechanism is reflected as platelet disappearance after they spend their 
designated transit time ( P L ) in the PL compartment. The second one is rather age- 
independent and it is reflected as constant platelet efflux (d) throughout all platelet 
age-stages. 
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TPO EFFECTS Treturn to Contents! 

Having thus formulated the compartment model for thrombopoiesis, we still 
have to define the effects of TPO on it (Fig.l). As was mentioned before we can 
divide the effects of TPO to three levels. These could be expressed in terms of 
dependence of the following parameters on TPO concentration: (i) amplification rate 
( )> (ii) the rate of cell maturation or, alternatively, transit time through given 
compartments ( ) , and (iii) the fraction of megakaryocytes of a given ploidy that 
undergoes an additional endomitosis and passes to the next ploidy class (). N e i thor 
For simpl icity we assume that neither the rate of platelet release per given 
megakaryocyte nor the rate of platelet utilization (d) is dependent on TPO 
concentration in our model . 3 * [RRFFTlFNPn] F W we shall define the 
pharmacokinetics of TPO. 



TPO concentration, [return to Contents] 

Many works have shown that recombinant human full-length TPO, as well as its 
truncated form rHuMGDF, are fully active biologically. [REFERENCE] Therefore, in 
our model we simply add exogenously administered recombinant protein to 
endogenously produced one in order to calculate actual TPO concentration (Q. 

The rate of TPO production in the main TPO production sites, i.e. liver and 
kidney, was shown to be constant under thrombocytopenia or thrombocytosis of any 
severity. 2 ' 14 ' 15 However, the level of TPO mRNA in the sites like bone marrow and 
spleen, where it is produced at the significantly lower rate than in the liver and 
kidney, can change considerably. 2, 14 Despite the apparent importance of spleen and 
bone marrow concentrations as megakaryopoiesis mainly occurs at these sites, it was 
shown that there is no significant difference between TPO levels in the bone marrow 
and in the peripheral blood. 17 We therefore assume that the bone marrow and spleen 
production of TPO is negligible, and that endogenous TPO is produced at a constant 
rate (/?). 

The main mechanism that controls TPO concentration in the blood is receptor- 
mediated TPO uptake (w). 2 ' 5 ' 17 ' 18 ' 19 > 20 

Another mechanism, we assume, is non-specific TPO-receptor-independent 
clearance of TPO (/). We suggest that this mechanism is rather insignificant in normal 
state, when receptor-mediated TPO binding, endocytosis, and degradation (w) remove 
most of the TPO (it was shown experimentally in mice) 30 . But when the amount of 
TPO (O rises significantly above the ability of the receptor pool to remove it, the 
non-specific clearance (I) becomes important. 

Thus the formula for TPO levels C at time i+l is: 



Cm = C i +p + x i -u i -l i p 9 x i9 u i9 l t ^0 C,. >0 

where Q is TPO concentration at the current hour (i); C*/ + y is approximation of 
the TPO concentration at the next hour (detailed below); p is TPO produced per hour 
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endogenously; is the addition to TPO concentration due to exogenous TPO 
administration; w, is TPO removed from the blood by receptor-mediated binding; /,■ is 
TPO cleared from circulation by non-specific mechanisms. 

We assume that receptor-mediated TPO clearance depends on the total number 
of TPO receptors (n) and on the ability of each receptor to uptake TPO (a): 

^2) u i =n i -a n n a>0 

where represents the receptor pool and a is TPO-clearing ability of the 
receptors, i.e. amount of TPO that each receptor removes per hour. 

Both, megakaryocyte and platelet mass contribute to the total receptor number 
(«) and, thus, to the rate of TPO clearance (w). 15 We assume that every platelet bears 
the same number of TPO receptors (mpi). The receptor number on megakaryocytes, 
however, changes. Thus, the receptor pool (n) is: 



\ ) ™ — . . - m 

comp, j, i comp 




+ IPL i * m PL <1 comp J ,i > 9pl 4 > ™ CO mpJ >™PL^ 0 



J 



where comp (1 to 4) is one of the platelet releasing megakaryocyte 
compartments (MK16, MK32, MK64, MK128, respectively) ;y is the period (in hours) 
that a given megakaryocyte already spent in the specific compartment; [] denotes 
rounded to an integer; q C om P j,i is the quantity of the megakaryocytes of the specific 
compartment {comp), which spent a given period (j) in it; m compJt is the receptor 
number on given megakaryocyte; q PL , is the platelet number; m PL is the receptor 
number per platelet. 

It is assumed that the number of TPO receptors on each megakaryocyte (m compJ ) 
equals the number of platelets that the megakaryocyte is capable of releasing (c comp ) 
times the average number of receptors per every potential platelet (b). 

m =(r —r •/)•/> r r i h>Ci 

'"compj yrcomp 'comp J ) comp ' comp ? 7 ? ^ — w 



where c comp is the number of platelets that the megakaryocyte of the specific 
compartment comp can release during its entire life-span ( co mp )\ r comp is the rate of 
platelet release by the megakaryocyte; j is the period that this megakaryocyte already 
spent in this compartment; b is the number of receptor on the megakaryocyte per 
potential platelet. 

We assume that the non-specific TPO clearance (//) is exponential, i.e. every 
hour some fraction (J) of a current amount of TPO (Q is removed from circulation: 



(5) 



/ f =/-C, f>0 



where / is the coefficient of non-specific TPO clearance and Q is the current 
TPO concentration. 
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Exogenous TPO is included in the model as a linear relation of the initial 
maximum TPO blood concentration (xi) to the administered intravenous (IV) dose (s) 
(the relation coefficient is 0.0167) 21 : 



(6) 



jc ; =0.0167-5, 



j, >0 



The state when TPO completely disappears from the blood seems very unlikely 
based on biological logic, so we restricted the lower limit of possible TPO 
concentration to certain minimum (positive). Thus, the equation (1) is modified to 
receive the full TPO concentration equation: 

(7) C M = max((C. +p + x i -a-n i -f-C i ) 9 e) £>0 
In steady state, the TPO concentration (C) is constant. 
TPO effects on amplification rate, [return to Contents] 

In our model, there are only two compartments, SC and CFU-Meg, whose cells 
are capable of dividing. These compartments differ significantly from each other, 
thus, we shall discuss them separately. Cells of other model compartments do not 
proliferate, and so their amplification rate () equals 1 under all circumstances. 

SC compartment, [return to Contents ! 

Since TPO is primarily a thrombopoiesis-stimulating cytokine, we assume that 
the cells, which are not committed to thrombopoietic line yet (the SC compartment in 
our model), are relatively insensitive to TPO, compared to committed megakaryocytic 
cells. In our model this is considered as a threshold () i n TPO concentration. Only 
above this threshold ( ) TPO affects stem cells. As long as TPO remains below the 
threshold ( ), stem cells in the model are regulated by intrinsic TPO-independent 
mechanism that keeps the size of their population almost constant. 

Thus, below the threshold () , SC amplification rate (sc) i s determined hourly 
depending on the current number of cells in the SC compartment. It seems 
biologically reasonable that the dependence equation is a sigmoidal function where 

sc changes from 1 (i.e., no amplification, the cell number remains the same) when 
the cell number approaches infinity, up to the maximal value sc w when the cell 
number approaches zero. The increase in amplification rate ( S c) i s relatively gradual 
as long as the cell number (q sc ,) exceeds certain critical value (we assumed it to be a 
fraction (v) of the normal cell number (gscnorm))- However, when the cell number falls 
bellow this threshold, S c begins to increase rapidly in order to restore the SC 
compartment as soon as possible. We assume that at normal cell numbers (qscnorm), 

s c i should be a fraction (y) of its maximal value S c w . Following is an example of 
such equation: 




norm 



a~sc=(a SCw -l)- 



f 



i \ 

1 v Sl +1 



1 



+ 1 



q sc ,S h2 ,0>O 
0<y,v<\ 



where * sc, /+/ is the amplification rate calculated based solely on the cell 
number; S c w is the maximal possible rate of cell amplification in the SC 
compartment when TPO concentration (C,-) is below the threshold (); q S c, i is a 
quantity of cells in the SC compartment; q S c norm is the normal quantity of cells there. 
Si and S2 are the sensitivity coefficients in the regions of qsa higher or lower than the 
critical value {vq$c norm), respectively. These values determine the sensitivity of the 
mechanism that links the amplification rate (s c) with the cell number (qsa)- In other 
words, they determine the steepness of the dependence curve in the corresponding 
regions. High Sj or S2 mean that $ c changes significantly due to slight changes of 
qsa and low Sj or S2 mean that 5 c remains relatively constant whatever the changes 
of qsc are. Distinguishing between Sj and S? allows us to force the amplification rate 
( sc) to grow rapidly as the cell number (qsa) falls below the critical value, thereby 
increasing the resistance of the system to further cell number (qsa) decay. Although 
the symbols Sj and S2 appear in several equations, their values are specific for every 
equation. 

We suggest that TPO concentration (C) increase above the threshold () should 
occur in severe platelet and/or megakaryocyte deficiency, or when TPO was 
administered exogenously. We assume that at these circumstances, TPO further 
increases the rate of cell amplification in the "Stem Cell" compartment ( S c)- W e 
assume also that the increase is proportional to the difference between actual TPO 
concentration (Q) and the threshold () . Thus, TPO effects appear gradually from the 
zero increase, when TPO concentration (Q) equals the threshold (). Saturation of 
the mechanisms of TPO effect is reflected in the concavity of the effect function. 

The following is an example of such a function: 

(9) 

<*sc, m {<lsc, i >Ci ) = a sc. w [Qsc, i > C s )+ 1 • ln(Q - 9 + 1) t 9 0 > 0 

C- t >0 



where S c, 1+1 is the same expression as in equation (8), i.e. amplification 
calculated based on a TPO-independent mechanism, and the second operand is the 
TPO-related contribution to the amplification rate (sc , /+/)• t determines the steepness 
of the dependence curve (t is non-negative). Although the symbol t appears in several 
equations, its value is specific for every equation. We add one to the In argument in 
order to ensure positivity of the In result. 

CFU-Meg compartment, rreturn to Contents ] 

In contrast to the cells of the SC compartment, we assume that cells of this 
compartment are fully sensitive to TPO and respond to the absolute TPO 
concentration (C,), not to its difference with a threshold (Q - ). In our model, they 
have no TPO-independent proliferative mechanism, and cease to proliferate when 
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deprived of TPO. On the other hand, when TPO concentration (Ci) in the system 
increases, cfu does not rise to infinity, but rather gradually reaches saturation, which 
also seems reasonable biologically. At normal TPO concentrations (C nor m), we assume 
cfu to be a fraction (h) of its maximal value (c Ft/max). Thus an equation that 
describes the relation of the amplification rate of CFU-Meg cells ( C f u) to TPO 
concentration (Ci) represents a sigmoidal function with C f u equaling 1 when TPO 
concentration (Ci) is zero, passing through h times cf c/max when TPO concentration 
is normal (Ci), and approaching an asymptote in cf t/max when TPO concentration 
(Ci) approaches infinity. In addition, in order to enable the system to be sensitive both 
to the regulation by endogenously produced TPO and to the effect of the exogenously 
administered drug, it was assumed that the function changes relatively rapidly in the 
region of normal TPO concentration (C n0 rm) and with a much smaller rate when a 
TPO concentration (Ci) is somewhat higher than normal (C norm ). The following is an 
example of such a function: 



(10) 
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where C fu, /+; is an amplification rate of the CFU compartment; cf Umax is a 
maximal value of amplification rate there; C nor m is normal TPO concentration; t is the 
parameter that determines the steepness of the dependence curve. 

TPO effects on transit time, [return to Contents] 

For the reason noted earlier, we assume that all platelet-releasing 
megakaryocyte compartments have the same transit time (mk). We assume that 
neither the relation of megakaryocyte volume (and thus, its platelet releasing capacity 
Ccomp) nor of its rate of platelet release r comp to megakaryocyte ploidy, is affected by 
TPO. Therefore, the transit time ( M k) through the noted compartments is constant. 
Platelets also spend in average a constant time in the circulation ( PL ), which is not 
affected by TPO concentration (Ci). 

In contrast, the transit times of the SC, CFU-Meg, and MKB compartments 
( sc, cfu, mkb, respectively) are functions of the micro-environmental conditions. 
Since cells that are far from maturation are not expected biologically to undergo a 
sudden shift to maturation, it seems that these functions determine the value the transit 
time should approach, rather than the actual transit time. The actual transit time 
changes gradually: every hour it changes by l-2eee hours towards the function- | 
determined value. Thus, we will mean the value, that transit time approaches rather 
than the transit time itself when speaking about transit time () c alculations below. 

SC compartment, [return to Contents ] 

We assume that regarding transit time ( s c ), the SC compartment differs from 
others in the same way as regarding amplification rate ($c)- It means that the cells of 
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this compartment respond to TPO only when its concentration (Q) rises above a 
threshold ( ). This threshold () is the same as for the amplification rate ( S c)- 
Below the threshold ( ) SC transit time ( s c) is assumed to be regulated by a TPO- 
unrelated mechanism dependent on the current cell number {q sc ,) only. The function 
of this dependence changes the transit time ( s c) from its minimal value ( sc u ) when 
the cell numbers (qsa) approach infinity, through the normal value that is greater than 
the minimal one by factor g, up to the highefst valuej Csc maxl * limited determined 
solely by biological reasons. This means that when the cell number in SC 
compartment (q S c ,•) is high, the cells will pass rapidly to the next compartment, thus 
reducing the SC one; and they will remain longer in the SC compartment when their 
number (qsc d is low, thus repopulating it. This manner of regulation seems 
reasonable biologically. 

We suggest that similar to the amplification rate (sc), the transit time ( S c) i n 
the range of very low cell numbers (q S c /) (lower than a certain fraction (v) of the 
normal (qscnorm)), is very sensitive to further cell number decrease, and grows rapidly, 
thereby resisting compartment exhaustion. This fraction (v) is the same as for the SC 
amplification rate. 

Following is an example of such a function: 
(11) 



TscMx{q SCJ ,c) = 



r 



SCu 



HsCnorm 



Qscj 



J 



> Qscj^ v 'Qsc 



norm 



J 



-(■ 



)■ 



f 



QsCJ 



QSC norm J 



> qscj <v '<isc 



norm 



T SC — T 



SCu 



V 



s h2 > 0 



T SCu >° 

0<v<l 



where sc, i+i is the transit time calculated based on cell numbers (qscj) ° n ly> 
i.e. when TPO concentration (Q) remains below the threshold (); S c u is the 
minimal possible transit time through SC in these circumstances; v is the fraction of 
normal cell number (qscnorm)> below which the dependence of the transit time Csc, 
on the cell number (qscj) changes; Sj and S 2 are the sensitivity coefficients in the 
regions of qscj lower and higher than vqscnorm, respectively. 

If TPO concentration (C/) in the model rises above the threshold (), the transit 
time ( sc) is assumed to shorten in a dose dependent manner. As for the amplification 
rate ( S c)> its decrease is presumed to be proportional to the difference between actual 
TPO concentration (Q and the threshold () . However, a shortening of the transit 
time ( sc) down to zero by TPO is biologically illogical , so we suggest that the 
transit time ( sc) approaches some minimal value as TPO concentration (C,*) 
increases. In our model this minimum represents a fraction (k) of the transit time 
calculated on the basis of cell numbers (*s c, /+/) as described earlier (equation (1 1)). 
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Following is an example of such an equation: 
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where 5C , /+/ is the transit time when TPO concentration (C,-) is higher than the 
threshold ( ); sc, i+i is the transit time calculated on the basis of cell numbers as 
described in equation (1 1); k is the fraction of * s c, i+i that gives the minimum, which 
the transit time approaches as TPO concentration (Q) increases; C determines the 
point of TPO concentration (C,), around which the transit time ( S c) i s the most 
sensitive to concentration (Q change; t determines the steepness of the dependence 
curve (t is non-negative). Multiplication by t enables td regulate the sensitivity to Q 
with t < 1 in the same manner as when t > 1 . 

CFU-Meg and MKB compartments, (return to Contents ] 

We assume that the transit time parameters of these two compartments ( C fu, 
mkb, respectively) are dependent solely on TPO and respond to the absolute TPO 
concentration (Q), rather than to its difference with a threshold (Q - ). As TPO 
level (C,-) drops, the cell passage through these compartments slows, i.e. transit time 
( comp) increases up to the values limited solely by biological reasons ( co m p tma x) (we 
suggest that the cells cannot stay in these compartments for an infinite period of time). 
On the other hand, when the TPO concentration (Q) in the system increases, com P 
does not shorten to zero, but rather asymptotically reaches co mPrm i n9 thus bounding the 
function from below. This also seems biologically reasonable, as the cells cannot 
move through the compartment in one instant. At normal TPO concentrations (C norm ), 
we set comp to equal its normal value ( c 0 mPt norm)- In addition, in order to enable the 
system to be sensitive both to the regulation by endogenously produced TPO and to 
the effect of the exogenously administered drug, it was assumed that the function 
changes relatively rapidly in the region of TPO concentrations (Q) lower than normal 
(Cnorm) and with a smaller rate when a TPO concentration (Q) is higher than normal 
(Cnorm)- The following is an example of such a function: 
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where comp is one of the aforementioned compartments (CFU-Meg or MKB); 
compj+i represents the transit times through these compartments; co mPtmin , com p.norm, 
and comp.max are the minimal, normal and maximal transit times when TPO 
concentration is normal; t\ and ti determine the steepness of the dependence curve in 
the regions of Q lower and higher than C norm , respectively. 

TPO effects on the fraction of cells that flow from one compartment to the 

next, [return to Contents] 

The discussed parameter is the proportion of "mature" cells that passes to the 
next compartment at any given moment ( ) . We designate as "mature" the cells that 
are potentially ready to pass to the next compartment but do not necessarily pass. As 
was noted earlier, we assume that the fraction of the SC compartment that commits to 
the megakaryocytic lineage ( s c) is constant and TPO-independent. From the two 
subsequent compartments, CFU-Meg and MKB, every mature cell is supposed to 
emerge to the next compartment in any case. Thus, TPO in our model does not affect 
in the first three compartments. 

In contrast, the fraction of MK16-, MK32-, and MK64- megakaryocytes that 
continues with additional endomitoses and flows to the next compartment ( com p ) is 
assumed to be in the range of 0 to 1 depending on TPO concentration (C,). Because 
there is no compartment with ploidy greater than 128N, the megakaryocytes of the 
MK128 compartment do not flow to any other compartment. 

The dependence of MK16, MK32, and MK64 parameters on TPO 
concentration assumed to be delayed with calculated based on TPO concentration 
that was before the period needed for one additional endomitosis (C,-_). The reason 
for this is the assumption that the cells that enter a given megakaryocyte compartment 
first "decide" what they do, i.e. undergo another endomitosis and not participate in 
platelet release in this compartment, or alternatively, begin with platelet release and 
remain in this compartment until complete exhaustion. Thus, the cells that leave a 
compartment now were preprogrammed before they start an additional endomitosis. 

In the model the dependence is expressed by a sigmoidal function with s e t to 
0 when TPO concentration (Q) is 0, equaling the normal value ( no rm ) when TPO 
concentration is normal (C norm ), and approaching 1 asymptotically. 
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where comp is one of the discussed compartments (MK16, MK32, or MK64); 
comp,i+i is a parameter of these compartments; i s the time needed for one 
additional endomitosis; c 0 m P ,norm is the value of co mp under normal TPO 
concentration (C norm ); t determines the steepness of the dependence function (t is non- 
negative). 
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The time needed for an additional endomitosis ()a s sumed to be the same in the 
three relevant compartments (MK16, MK32, and MK64). 

PROGRAM STRUCTURE [return to Contents! 

Our model was built while trying to imitate real bone marrow physiology. Each 
compartment is subdivided into small sections that contain the cells of a specific age 
with a resolution of one hour. For example, the fifth age-section of the MKB 
compartment contains the cells which spent 5 hours in this compartment. Every hour, 
all the cells in the virtual bone marrow pass to the next age-section in the same 
compartment. 

When the cells spend all the transit time predetermined for them in a given 
compartment, they pass to the beginning of the next compartment, i.e., every hour the 
cells that have just left each compartment fill the first age-section of the next one. The 
cells that leave the MK128 compartment disappear. The first age-section of the SC 
compartment is given a constant fraction of the number of cells exiting the 
compartment. The platelet releasing cells contribute platelets to the first age section of 
the PL compartment every hour. 

The described structure is implemented as an array of 8 rows and 360 columns 
(Fig. 2). The 8 rows relate to 8 cell compartments, and the columns to the age 
sections, as we assumed that none of the transit time could exceed 360 hours. This 
array which will later be referred to as age-distribution, is refilled hourly according to 
the rules described above. The smoothened graphic representation of this array is 
shown in Fig. 3. Although the age-distribution array is useful in model investigation, 
this is not what physicians and biologists are interesting in. Rather the time course of 
the sizes of different compartments is required. Hence, in our simulation we file the 
hourly states of all the compartments and TPO concentration at all times. For this 
purpose, another array, which will later be referred to as time course is established in 
the model (Fig. 4). The rows in the array are cell compartments and the columns are 
the time steps of the simulation course. Every time-step of the simulation (one hour of 
"patient's life"), the number of cells in all age-sections is summarized for each 
compartment and the next column in the time-course array is filled. Thus, every cell 
in this array represents the total number of cells in a given compartment at a given 
time point. 

There is an additional row in the time-course array that relates to the TPO 
concentration in the blood. TPO concentration is written down every time-step 
concurrently with the cell numbers. Example of it is shown in Fig. 5. 

Parameter determination [return to Contents ] 

Our model is very flexible allowing the user to fit it to patients with diverse 
blood and bone marrow parameters. Real patients differ in their baseline platelet 
counts (gpinorm)~ and numbers of bone marrow precursors (q C ompnorm), in the | 
sensitivity of their stem cell "intrinsic" regulation mechanism (S-s), in their minimum 
and normal transit times (s c,u, c omp.norm) and maximal amplification rates (sc, w, 
comp.max), rates of platelet release by megakaryocytes (r comp ), fractions that each 
megakaryocyte ploidy class contribute for additional endocytosis ( co mp.norm), and the 
time needed for endomitosis 0 . Furthermore, the baseline TPO level (C norm ), the rate 
of TPO production (p), receptor- and non-receptor-mediated TPO clearance (a and f), 



15 



the threshold of TPO effect on the SC compartment (), and the sensitivity of 
different cell parameters to TPO also differ between patients (/-s).— 

Since the aim of this work is to achieve the best fit of the model to each patient, 
the patient-related parameters should be determined individually for each patient. 
However, practically, it would be very difficult to predetermine many of these 
parameters for every patient. We calculated certain average parameters (Table 1), 
based, generally, on published data (see references in the Table). These averaged 
parameters are used as a framework, into which known individual characteristics are 
included. Thus, before the simulation of the specific patient is begun, the user should 
correct the parameters of the model as far as he knows the relevant information about 
the patient. 

Usually, the known patient-related data are not parameters in the form defined 
by our model, but rather some kind of measurements obtained in the clinic (e.g., day 
and value of post-chemotherapy thrombocytopenia nadir, day and value of platelet 
peak after TPO administration, or even the change in the megakaryocyte modal ploidy 
following some perturbation). [REFERENCE] In this case, we analyze the available 
data to convert it into the model-compatible format. 

Sometimes, the only available patient-related data is the graphic representation 
of patient's platelet course following some perturbation (e.g., cell-suppressive therapy 
or TPO administration). [REFERENCE] The data may also be a picture of the platelet 
course without any external perturbation (e.g., cyclic thrombocytopenia). 
[REFERENCE] In these cases we change the model parameters by trial-and-error 
until a good compliance of the model graphic output with the given patient's graphs is 
achieved. However, some kind of analysis is possible even in these cases, and the 
choice of the parameter sets for trying is not absolutely random. 

Results [return to Contents] 

We implemented the described model as a computer simulation that retrieves 
the changes that occur in the human thrombopoietic system (platelet counts, bone 
marrow precursor numbers, and TPO concentrations) over several years. The 
resolution of the simulation output is one hour. Time units and periods that we will 
speak of relate to the simulated patient's life, not to the running time of the program. 
The programming tools give the user the following features: 

1) As was noted in the previous section, the user can set the baseline values 
and all other thrombopoietic parameters as far as he knows them about the 
patient before running the simulation. 

2) The user (e.g., physician, scientist) can determine how long the simulation 
will run, from 12 hours up to several years. 

3) The user can determine how frequently he wants the program to show him 
the course of a patient counts up to the moment. The frequency can change 
from every 12 hours to once during the overall period of the simulation. 

4) The user can determine what resolution he would like to see the output 
graph, one hour or above. 
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5) The user can determine if he wants to see the graphical representation of the 
age distribution through the compartments at any moment of the simulation. 

6) The user can simulate a cell-suppressive therapy at any moment while 
running the simulation by reducing one or several of the compartments by 
any value. 

7) The user can simulate exogenous TPO administration at any moment while 
running the simulation by controlling: 

• dose height 

• number of dosings 

• frequency of dosings 

The simulation tool has been carefully tested with respect to the published 
experimental results, and has proved to be well calibrated for averaged human 
thrombopoiesis. The technique for relatively fast modifications in order to fit a 
specific patient has been also developed. Fig. 6 A and 7 A show examples of the tool's 
output. Here we retrieve the in vivo platelet counts following TPO administration 
(Fig. 6A) and chemotherapy without TPO (Fig. 7A). 

Discussion [return to Contents] 

We built a model that simulates cell and platelet counts in the steady state, as 
well as after perturbations to the hematopoietic system, e.g., cell-suppressive therapy, 
recombinant thrombopoietm administration, etc. It is possible to simulate any protocol 
of drug administration and any hematological state of a patient, regarding his/her 
platelet count and number of bone marrow megakaryocytes and their precursors. The 
model is very flexible and can be adapted to many categories of patients, or healthy 
platelet donors. It can also be modified to fit species other than human. By providing 
specific parameters one can adjust the model so as to yield particular predictions 
about the thrombopoietic profile of an individual patient. We believe that platelet 
disorders, such as cyclic thrombocytopenia, can also be simulated in our model. 

This model makes it possible to check any given treatment protocol and to 
choose the best one according to the user's criteria. As such, it strongly reduces the 
number of necessary clinical trials and lowers the risk of treating a patient with 
ineffective TPO protocol. 

In the near future, we will operate an optimization program on our model. For 
any individual this modified program will select the best protocol out of an entire 
spectrum of possibilities. The criteria by which the optimization module determines 
an optimal treatment protocol may be one or several of the following (this list is just 
an approximation): 

1 . The smallest number of TPO dosings required for achieving any given 

aim. 

2. * The lowest maximal TPO dose required for achieving any given aim. 
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3. The minimal total amount of TPO needed for rehabilitation. 

4. The smallest deviation from the baseline at platelet nadir after 
chemotherapy or another cell-suppressive treatment. 

5. The shortest period of thrombocytopenia. 

Fig. 6 and 7 show examples of the tool's output. Here we retrieve the 
experimentally obtained in vivo platelet counts following TPO administration (Fig. 
6A) and chemotherapy without TPO (Fig. 7A). Most notably, we show that by use of 
our tool, it is possible (at least theoretically) to construct a TPO protocol that would 
reduce by 75% the amount of TPO currently used for chemotherapy patients, while 
maintaining the same thrombocyte rehabilitation pattern. 



References [return to Contents ] 

1 . Harker LA, Hunt P, Marzec UM, Kelly AB, Tomer A, Hanson SR, Stead RB: Regulation of 
platelet production and function by megakaryocyte growth and development factor in nonhuman 
primates. Blood. Mar 1996; Vol. 87(5); pp. 1833-1844. 

2. Alexander WS: Thrombopoietin. Growth Factors. 1999; 17(1); pp. 13-24. 

3. Kaushansky K: Thrombopoietin: the primary regulator of platelet production. Blood. Jul 1995* 
Vol. 86(2); pp. 419-431. 

4. Vadhan_Raj-Raj S: Recombinant human thrombopoietin: clinical experience and in vivo biology. 
Seminars Hem. Jul 1998; Vol. 35(3); pp. 261-268. 

5. Harker LA: Physiology and clinical applications of platelet growth factors. Current Opinion 
Hematol 1999; Vol. 6; pp. 127-134. 

6. Levin J: Thrombopoietin - Clinically realized? New Eng J Med. Feb 1997; Vol. 336(6); pp. 434- 
436. 

7. Somlo G, Sniecinski I, ter Veer A, Longmate J, Knutson G, Vuk-Pavlovic S, Bhatia R, Chow W, 
Leong L, Morgan R, Margolin K, Raschko J, Shibata S, Tetef M, Yen Y, Forman S, Jones D, 
Ashby M, Fyfe G, Hellmann S, Doroshow JH: Recombinant Human thrombopoietin in 
combination with granulocyte colony-stimulating factor enhances mobilization of peripheral blood 
progenitor cells, increases peripheral blood platelet concentration, and accelerates hematopoietic 
recovery following high-dose chemotherapy. Blood. May 1999; Vol. 93(9); pp. 2798-2806. 

8. Neelis KJ, Hartong SC, Egeland T } Thomas GR, Eaton DL, Wagemaker G: The efficacy of single- 
dose ao^runistration of thrombopoietin with coadministration of either granulocyte/macrophage or 
granulocyte colony-stimulating factor in myelosuppressed rhesus monkeys. Blood. Oct 1997* Vol 
90(7); pp. 2565-2573. 

9. Murray LJ, Luens KM, Estrada MF, Bruno E, Hoffman R, Cohen RL, Ashby MA, Vadhan-Raj S: 
Thrombopoietin mobilizes CD34 + cell subsets into peripheral blood and expand multilineage 
progenitors in bone marrow of cancer patients with normal hematopoiesis. Exp Hem. 1998- Vol 
26; pp. 207-216. 

10. Basser RL, Rasko JE, Clarke K, Cebon J, Green MD, Grigg AP, Zalcberg J, Cohen B, O'Byrne J, 
Menchaca DM, Fox RM, Begley CG: Randomized, blinded, placebo-controlled phase I trial of 
pegylated recombinant human megakaryocyte growth and development factor with filgastrim after 
dose-intensive chemotherapy in patients with advanced cancer. Blood. May 1997; Vol 89(9)* pp 
3118-3128. 

11. Gordon AS : Regulation of hematopoiesis. N.-Y. 1970 Vol. 2, Section IX (textbook). 

12. Harker LA, Finch CA: Thrombokinetics in man. J Clin Invest. 1969; Vol.48; pp. 963-974. 



18 



• 



13. Harker LA: Thrombokinetics in idiopathic thrombocytopenic purpura. Br J Haematol. 1970; Vol. 
19; pp. 95-104. 

14. Sungaran R, Markovic B, Chong BH: Localization and regulation of thrombopoietin mRNA 
expression in human kidney, liver, bone marrow, and spleen using in situ hybridization. Blood. Jan 
1997; Vol. 89(1); pp. 101-107. 

15. Nagata Y, Shozaki Y, Nagahisa H, Nagasawa T, Abe T, Todokoro K: Serum thrombopoietin level 
is not regulated by transcription but by the total counts of both megakaryocytes and platelets 
during thrombocytopenia and thrombocytosis. Thromb Haemost. 1997; Vol. 77; pp. 808 - 814. 

16. Nagahisa H, Nagata Y, Ohnuki T, Osada M, Nagasawa T, Abe T, Todokoro K: Bone marrow 
stromal cells produce thrombopoietin and stimulate megakaryocyte growth and maturation but 
suppress proplatelet formation. Blood. Feb 1996; Vol. 87(4); pp. 1309-1316. 

17. Hsu HC, Tsai WH, Jiang ML, Ho CH, Hsu ML, Ho CK, Wang SY: Circulating levels of 
thrombopoietic and inflammatory cytokines in patients with clonal and reactive thrombocytosis. J 
Lab Clin Med. 1999; Vol. 134(4); pp. 392-397. 

18. Stoffel R, Wiestner A, Skoda RC: Thrombopoietin in thrombocytopenic mice: evidence against 
regulation at the mRNA level and for a direct regulation role of platelets. Blood. Jan 1996; Vol. 
87(2); pp. 567-573. 

19. Alexander WS: Thrombopoietin and the c-Mpl receptor: insights from gene targeting. IntJ 
Biochem Cell Biol 1999 Oct; Vol. 31(10); pp. 1027-1035. [ABSTRACT] 

20. Miyazaki M, Fujiwara Y, Isobe T, Yamakido M, Kato T, Miyazaki H: The relationship between 
carboplatin AUC and serum thrombopoietin kinetics in patients with lung cancer. Anticancer 
Research. 1999; Vol. 19; pp. 667-670. 

21. Vadhan-Raj S, Murray LJ, Bueso-Ramos C, Patel S, Reddy SP, Hoots WK, Johnston T, 
Papadopolous NE, Hittelman WN, Johnston DA, Yang TA, Paton VE, Cohen RL, Hellmann SD, 
Benjamin RS, Broxmeyer HE: Stimulation of megakaryocyte and platelet production by a single 
dose of recombinant human thrombopoietin in patients with cancer. Ann Intern Med. May 1997; 
Vol 126(9); pp. 673-681. 

22. Wichmann HE, Gerhardts MD, Spechtmeyer H, Gross R: A mathematical model of 
thrombopoiesis in rats. Cell Tissue KineU 1979; Vol. 12; pp. 551-567. 

23. Harker LA, Roskos LK, Marzec UM, Carter RA, Cherry JK, Sundell B, Cheung EN, Terry D, 
Sheridan W: Effects of megakaryocyte growth and development factor on platelet production, 
platelet life span, and platelet function in healthy human volunteers. Blood. 2000 Apr; Vol. 95(8); 
pp. 2514-2522. 

24. Mayani H, Dragowska W, Lansdorp PM: Lineage commitment in human hemopoiesis involves 
asymmetric cell division of multipotent progenitors and does not appear to be influenced by 
cytokines. J Cellular Physiol 1993; Vol. 157; pp. 579-580. 

25. Golde DW: The Stem Cell. Medicine. Dec 1991. 

26. Morrison SJ, Uchida N, Weissman IL: The biology of hematopoietic stem cells. Annu Rev Cell 
DevBiol 1995; Vol. 11; pp. 35-71. 

27. Eller J, Gyori I et ah Modelling thrombopoiesis regulation - 1: model description and simulation 
results. ComputMath Applic. 1987; Vol. 14 (9-12); pp. 841-848. 

28. von Schulthess GK, Gessner U : Oscillating platelet counts in healthy individuals: experimental 
investigation and quantitative evaluation of thrombopoietic feedback control. Scand J Haematol. 
1986; Vol. 36; pp. 473-479. 

29. Harker LA, Marzec UM, Hunt P, Kelly AB, Tomer A, Cheung E, Hanson SR, Stead RB: Dose- 
response effects of pegylated human megakaryocyte growth and development factor of platelet 
production and function in nonhuman primates. Blood. Jul 1996; Vol. 88(2); pp. 511-521. 

30. Fielder PJ, Gurney AL, Stefanich E, Marian M, Moore MW, Carver-Moore K, de Sauvage FJ: 
Regulation of thrombopoietin levels by c-mpl-mediated binding to platelets. Blood. Mar 1996; 
Vol. 87(6); pp. 2154-2161. 



19 



31. Swinburne JL, Mackey MC: Cyclical thrombocytopenia: characterization by spectral analysis and 
a review. J Theor Medicine. 2000; Vol. 2; pp. 81-91. 

32. Rasko JEJ, Beglev CG: Molecules in focus: The thrombopoietic factor. Mpl-ligand. IntJ Bioch 
CellBioL 1998: Vol. 30: 657-660. 

33. De Sauvage FJ. Carver-Moore K. Luoh S-M. Ryan A. Dowd M. Eaton PL, Moore MW: 
Physiological regulation of early and late stages of megakaryocytopoiesis bv thromboooietin. J 
Esp Med. 1996 Feb: Vol. 183: 651-656. 

34. Tacke F. Schoffski P. Trautwein C. Martin MIL Stangel W. Seifried E. Manns MP. Ganser A. 
Petersen D: Endogenous serum levels of thrombopoietic cytokines in healthy whole-blood and 
platelet donors: implication for plateletpheresis. Br J Haematol. 1999: Vol. 105: 511-513. 

35. Beutler E. Lichtman MA. Coller BS. Kipps TJ: Williams HEMATOLOGY. 5 th edition McGraw- 
Hill. Inc. 1995: Chapter 118: pp. 1149-1161. 

36. Schofield R. Lord BI et ah Self-maintenance capacity of CFU-S. J Cellular Phisiol. 1980: Vol. 
103: 355-362. 

37. Rosendaal M. Hodgson GS. Bradley TR: Organization of haemopoietic stem cells: the generation- 
age hypothesis. Cell Tissue Kinetics. 1979: Vol. 12: 17-29. 



20 



APPENDIX B: 

A Model of the Neutrophil Bone Marrow and Peripheral Blood Compartment under the 
Effects of Growth-Factors and its Use as a Tool for Optimizing Treatment with 

Granulocyte Colony Stimulating Factor (G-CSF) 



Abstract 

We present a mathematical model of the neutrophil compartment of the hemopoietic system. 
Our aim is to use this model as a basis for optimizing treatment with Granulocyte Colony 
Stimulating Factor (G-CSF) for patients undergoing myelosuppressive treatment. Therefore, 
we formulated the model as a piecewise linear system of equations, and implemented it in 
AMPL (A Mathematical Programming Language), thus guaranteeing the ability to use linear 
programming techniques in the treatment optimization effort. 

Underlying our model are two main sources of information and validation: Firstly, the 
mechanism implemented in the model is in accordancejw ith^sugg^tions and experimental 
evidence available in the literature regarding the kinetics of the lineage. However, these 
suggestions are often broad-brush strokes allowing a considerable degree of freedom 
especially in view of conflicting evidence. Therefore, we used as a second source of 
information the data available in the literature regarding the effects of G-CSF on blood 
n eutrophil lev els, and insisted that our predictions be in close accordance with these data. 
This ensures that our model is in close accordance with both the qualitative and the 
quantitative aspects of the experimental data. Both aspects are essential for correctly 
predicting optimal treatment protocols. 

Introduction 

The neutrophil lineage originates in pluripotent stem cells that proliferate and become 
committed to the neutrophil lineage. These cells then undergo gradual maturation 
accompanied with proliferation through morpho logically distinguishable mitotic 

comp^tments: Myeloblasts, promyelocytes and myelocytes. The myelocytes then mature 
and lose their capacity to proliferate, andlhuTenter the posFmitotic compartment. In the 
post-mitotic compartment the cells continue their gradual maturation, which is not 
accompanied with proliferation through the three morphologically distinguishable sub- 
compartments: Metamyelocyte, band and segmented-neutrophils. Cells exit the various sub- 
compartments in the post-mitotic compartment and enter the blood as neutrophils. They then 
migrate from the blood to the tissues. 



The Granulocyte-Colony Stimulating Factor (G-CSF) effects an increase in blood neutrophil 
levels primarily by increasing production in the mitotic compartment and shortening the 
transit time of the post-mitotic compartment. REVIEW G-CSF USAGE. . . 

Thus, the first compartment of the mitotic pool (myeloblast) receives an inflow of cells from 
st em-cell pr ecursors. The inflow of the other compartments receives the outflow from the 
previous one, subject to multiplication factors due to cell replication in the mitotic stages. For 
a review of granulopoiesis and experimental methods see Mary - 1984. 

Models regarding granulopoiesis in normal humans and in humans with pathologies of the 
bone marrow were suggested previously in order to give a coherent description of the kinetics 
of granulocytes from experimental data (Cartwright - 1964, Mary - 1978, Rubinow - 1974). 
In recent years Schmitz et al. developed a kinetic simulation model for the effects of G-CSF 
on granulopoiesis (Schmitz - 1993), and used it for the analysis of administration of G-CSF 
to patients suffering from cyclic neutropenia (Schmitz - 1995). However, the data Schmitz 
rests upon for his model has been more accurately assessed in recent years by Price et al. 
(1998) and Chatta et al. (1996). Actual empirical data regarding compartment sizes and their 
transit times was not incorporated into their model despite the importance of these data 
(Dancey et al. 1976). 

In this work we propose a mathematical pipeline model o f the neutrophil lineage with and 
without the effects of G-CSF which rests on the latest experimental data available. Our work 
is the first of its kind, as far as we know, to have formulated the mathematical model in a 
way that w ould facilitate the process of finding optimal protocol s for neutropenia induced by 
chemotherapy. We then formulate an objective func tion for reconstitution of a patient 's bone 
marrow, and implement the optimization model as a linear programming problem. This will 
allow us to use it as a tool for predicting the optimal protocol for G-CSF administration to 
patients undergoing myelosuppressive chemotherapy. 

In this article we shall present the mathematical model and its predictive power and discuss its 
possible applications for determining G-CSF treatment protocols for patients undergoing 
myelosuppressive treatment. 

Tools and Methods 

Since our long-term aim is to provide a tool for calculating opthnizedjreat^ we 
formulated the model with piecewise linea r^gi^doj^ndimplemented it in AMPL Plus 
software. AMPL is a comprehensive algebraic modeling language for problems in linear, 
nonlinear and integer programming. We used the AMPL language processor, which supports 
all the commands and syntax of AMPL. 




The Mathematical Model 
1. G-CSF 

The effects of G-CSF on the neutrophil lineage are relayed in the model in three stages. The 
first is the administered amount of cytokine given at time t, which is marked: G adm ^JU^ 1 

The Gadm vector serves as the control variable of the subsequently described optimization 
programme, as the main objective of this work is the optimization of G-CSF administrati on. 

The second stage represents the phaimacoki netic behaviour of G-CSF in cir culation. It 
incorporates for instance the half-life of G-CSF, and could in the future be modified to 
express more of the effects of time on G-CSF activity. This level is marked G* 

blood 

G-CSF is eliminated from the blood in a Poissonic manner (Stute, 1992): 

G blood =G 'blood ( 1 ~ _ ^ _ ) + G^dm (Eq. 1 . 1 ) 

1 Yl 



Where ty % is the half-life of G-CSF in the blood, and G x blood = G^ m 

We note, that recently data regarding the dependence of the half-life of G-CSF on the 
neutrophil counts has emerged (Terashi et al. 1999). In the absence of exact kinetics of G- 
CSF effects on the neutrophil lineage, we consider the half life to be constant, though this 
could easily be modified should more exact information emerge. 

Note, that we set as zero the effects and level of the endogenously produced G-CSF, and that 
only exogenously produced G-CSF is considered to effect the kinetic parameters. If more 
empirical data regarding the production of endogenous G-CSF is made available, we could 
obviously incorporate it in this equation. 

The third and final stage models the pharmacodynamic effects of G-CSF on the kinetic ^ 
parameters. As will be elaborated subsequently the dependence of the various kinetic 
parameters of the neutrophil lineage on the level of G-CSF in the blood is assumed to be 
through either non-decreasing concave or non-increasing convex functions. This reproduces 
the effects of faturation that are seen in clinical studies on the effects of G-CSF (Duhrsen 
1988), i.e. addition of G-CSF carries a lesser effect when its level in circulation is already 
high. 



Mitotic Compartment 

The mitotic compartment is divided into subcompartments. The kth subcompartment contains 
all cells of chronological age between k-l and k hours, relative to the time of entry into the 
mitotic compartment. The number of cells in subcompartment k at time t is marked as m[. 

k E {L.t} 

4=Z\(Gbiood) (Eq.1.2) 

where T is the transit time of the entire mitotic compartment, x is assumed to be the same and 
constant for all cells entering the mitotic compartment. 



Where h is a vector reflecting the flow of newly committed cells into the mitotic 
compartment. The biol ogical grounds for this definition is the existence of a myeloid stem 
cell reservoir, which is known to supply new committed cells to the mitotic compartment. 
However, the reservoir's actual kinetics are not very well explored empirically. We therefore 
fix /; to levels such that the overall size of the mitotic compartment as w ell as the kinetics of 
the neutrophils in circulation would match those obtained empirically. 

Obviously, any new biological data could help define these kinetics more accurately within 
the framework of this model. We note that our results indicate that the assumption of a 
constant rate of stem cells flowing into the mitotic compartment in the absence of G-CSF is 
plausible. 

For every n e {1. . d], and for every t y amplificatio n occurs at the exit from m l n : 



where: 

a n is a non-decre asing concave functionofG-CSF levels in the blood, which determines 
the factor of amplification in the hourly subcompartment n. If for instance no amplification 
occurs at subcompartment n 0 at time then 

a =1 

The size of the mo rphological sub-compartments in the mitotic compartment at time t is 
determined as: 

(Eq.1.4) 

n=n 0 

Where n 0 is the first hourly sub-compartment of a morphological sub-compartment and n } is 
its last hourly sub-compartment. Note, that the division into the morphological sub- 



compartments will be used only for fine-tuning of the kinetic parameters with the use of 
experimental data. 



The post mitotic compartment at time t is a single quantity of cells p\ such that: 



(Eq.1.5) 



h is a convex, non-increasing function of G-CSF levels in the blood, which takes values in 
the range of [0-1]. This definition entails p*>0. 

We shall mark as o t the outflow from the post mitotic compartment: 



The number of neutrophils in the circulating blood compartment at time t is marked ri and is 
modeled as a single quantity of cells, such that: 



where U/ 2 is the half-life of neutrophils in the blood, as defined in the biological literature. U/ 2 
is assumed to be held constant regardless of G-CSF levels (Lord 1989), though this could be 
easily modified. The kinetics of neutrophils in the tissues are not modeled in this work. 

This model will be incorporated into an optimization scheme that will have as its objective 
function both the aims of minimizing G-CSF administration and returning the neutrophil 
lineage to its normal levels. 

At G blood = 0> P = n, m T at ^ nor mal healthy level we have the following obvious 
relationship: 

y = "4 = = -^-xln2(Eq.l.8) 

Which reflects the stability of the steady state. 




(Eq.1.6) 



n' +1 =</+«Yl-J52 j 



(Eq.1.7) 



2. Linear Programming Implementation 

The above outlined model may be implemented in any number of optimization methods 
We chose to use linear programming because of its inherent advantages compared with other 
techniques, i.e. its ability to provide an o ptimal solution using at- least partially anal ytical 
methods, and therefore being more computationally tractable (Gill 1991). On the other hand, 
implementation of this model in linear programming carries with it the disadvantage that 
certain computations must be approximated linearly since they cannot be performed directly 
using linear methods. We shall therefore have to judge the 'closeness' of the solution 
obtained through linear programming with that obtained through another method of 
optimization that is not constrained with linear formulation. 

The significant points in the model that have to be modified due to the linear programming 
implementation are wherever multiplication of two variables is defined, since this operator is 
not itself linear. We therefore define multiplication as an approximated value constrained 
within piecewise linear constraints that most closely bound the product wit hin a four-fa ced 
polyhedron in 3- dimensional space whose vertices are 

{ (^miii > y min ' -*min * ^min )' (^min ' ^max ' *min * Xnax )> (^max ' ^min » ^max * ^min )' (^max > Jmax > **max * ^max ) 1 

Where x min , x max , y min > y max are the constant biologically defined minima and maxima of x and 
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(Eq. 2.1) 



Note, that multiplication may also be approximated with variations on the linear least squares 
method, by finding one plane that's closest to the four vertices defined (this should be Jested). 

The other functions that need to be defined linearly are those concerning the 
pharmacodynamics of G-CSF. Due to the nature of these functions (either non-increasing 
convex or non-decreasing concave), we implement these effects as piece-wis e linear functions 
whose breakpoints are the doses for-which actual experimental data is available (Chatta 
1994). Note, that the effects of G-CSF on each of the kinetic parameters were not determined 
in a detailed manner by experimentalists as far as we know. Rather its effects over a few dose 
levels on the neutrophil blood counts and the size of the morphologically different mitotic 
compartments and the post mitotic compartment was determined. From these we 
reconstructed the effects of G-CSF on the actual kinetic parameters (probability of mitosis, 
transit time and inflow of cells into the myeloblast compartment from sOtem cell progenitors) 



at the dose levels available in the literature, and then connected these points linearly to obtain 
piecewise linear functions relating G-CSF levels to their effect on those parameters. 
Obviously, more experimental data could be used to produce more accurate functions. 

We define a number of points in the mitotic compartment (see main text) as being potential, 
am plificatio njoints. At these points we use the linearly approximated multiplication operator 
(Eq 2.1) instead of the actual product defined in Eq. 1.3. 

At points where no amplification occurs the quantity from one compartment is simply 
transferred to the next: 

4tl =™n 

rrQ is set according to the steady state values of the mitotic compartment or is depleted 
according to the kill function of the chemotherapy. 

The flow out of the post mitotic compartment (Eq. 1.6) is similarly defined as a linear 
approximation of a product. 

3. Formulation of the Model as an Optimization Problem for Linear Programming 

The simulation spans a fini te number of d iscrete time steps denoted by T. 

We define as control variable the vector that represents G-CSF administration at every given 
hour ' : Gt adnt tE ^ 

The objective function is^efined as maximization of the fo llo wing expressio n: 

T 

^ (/?' ' P T ~ G X adm )(Eq.3.1) 
t = \ 

where p is the number of cells in the post mitotic compartment at time t, and • is a scalar 
weighting coefficient. The logic for formulating the objective function this way is that the 
ability to maintain the post mitotic compartment's steady state size for a prolonged period of 
time is sufficient for rehabilitation of the neutrophil lineage as a whole. Also, our goal is to 
minimize the total administered quantity of G-CSF. • is introduced to allow us to factor in 
both these goals in one objective. Also, this would allow a different weight to be given for 
certain times, e.g. were it determined (by clinicians) that the later states of the post-mitotic 
compartment should be weighted more than the first ones. Obviously this is only one of the 
possible formulations of the objective function as defined in the previous section. 



The pharmacokinetics and pharmacodynamics of G-CSF that were defined generally in the 
mathematical model are defined piecewise linearly. Some of the considerations that we put 
into formulating these functions were based directly on experimental evidence (elaborated in 
the main body of text). We note however, that actual experimental data regarding the direct 
effects of G-CSF on the kinetic parameters in which this model is interested is rather scant. 
Therefore, some formulations were conducted through partly analytic and partly trial-and- 
error methods. 



Determination of Parameter Values 
The Mitotic Compartment 

The mitotic compartment is where any long-ter m effects of G- CSF administration take place . 
Obviously, the post mitotic compartment's shortening of transit time due to G-CSF 
administration is the major contributor to heightened blood neutrophil counts in the short 
term. However, this high level cannot be maintained at all over the longer term without 
increased production in the mitotic compartment. 

The mitotic compartment was modeled with an intention to facilitate the specific cell-cycle 
cytotoxic effects of chemotherapy. Therefore, cohorts of one hour are modeled as undergoing 
a process of maturation and amplification (Eq. 1.3) culminating in their entry into the post- 
mitotic compartment (Eq. 1.5). Chemotherapy will be incorporated into the model by 
mapping the various cell-cycle phases (Gl, S, G2, M) to the hourly cohorts modeled and 
formulating a function of the cytotoxic effects of chemotherapy on these phases. 

The experimental literature shows wide agreement regarding the steady state normal sizes of 
the circulating neutrohpils, post-mitotic compartment and the three morphologically distinct 
sub-compartments of the mitotic compartment as well as post-mitotic transit time and the 
amplification rates in the mitotic sub-compartments (Dancey 1976, Price 1996, Dresch 1986). 
To determine other kinetic parameters of importance to our model, which were either not 
available in the literature or were given a wide range by experimentalists, we assumed steady 
state kinetics and employed an iterative process. These parameters include: the inflow of 
stem cells to the myeloblast compartment and the transit times of the mitotic sub- 
compartments. 

The half life of blood neutrophils and the steady state number of neutrophils was taken as 
7.6h and CUxlO 9 cells/kg body weight, respectively (Dancey, 1976). Similarly, the same 
calculation may be made for each patient that is to be modeled. This would allow the 
dynamics of every patient to be described by the simulation. We note, that the average size of 
the post-mitotic compartment (5.84xl0 9 cells/kg body weight - Dancey, 1976) and the transit 
time of the compartment (160h - Dancey, 1976; Mary, 1986; Price, 1996) are compatible 
with the size and half-life of the circulating neutrophil compartment reported by Dancey. This 
supports our steady state analysis. 

In order to determine the amount of cells in the hourly sub-compartments in the mitotic 
compartmen t we used a steady state assumption on all compartm ents in the lineage. The 
number of cells exiting the circulating neutrophil pool equals the number of cells exiting the 
post mitotic compartment (this was ascertained as explained above), which in turn equals the 



hourly production of cells in the mitotic compartment. This means that the number of cells in 
the last hourly cohort of the mitotic compartment can be determined from the neutrophil 
decay rate, which is available in the literature. We note however, that this calculation is based 
on the assumption, that there jsno ap optosis in the post-m itotic compartment. According to 
direct experimental data (Thiele 1997) it seems that apoptosis is not a significant phenomenon 
in the normal human bone marrow. If it is determined, that apoptosis occurs in the post 
mitotic compartment, our values for the production of cells in the mitotic compartment will 
obviously need to be modified upwards. We note however, that the size we calculate for the 
mitotic compartment is extremely close to that experimentally obtained by Dancey and Price, 
thus supporting the notion, that apoptosis is not a significant kinetic factor in the lineage. 

Regarding the transit time of the mitotic compartment there is little agreement in the literature 
with a range of 90-160 hours given by most experimentalists (Dresch, 1986). In order to 
determine the transit times of the mitotic morphological (Eq. 1.4) sub-compartments we 
considered the following constraints: The sizes of the theoretically obtained morphological 
sub-compartments must fit those reported experimentally in normal human hematopoiesis 
(Dancey, 1976) and under the effects of G-CSF (Price, 1996); at least 24h -the typical cell 
cycle- must separate amplification points; the size of the last hourly sub-compartment must 
equal the hourly production of the mitotic compartment (calculated with the aforementioned 
iterative process assuming steady state kinetics); amplification inside the compartment is set 
at the levels determined by Mary 1984; and finally, the total transit time of the mitotic 
compartment must be within the 90-160h range. We obtained an excellent fit with all of these 
constraints by using the values given in Table 1. 

Note, that when attempting other alternatives with shorter transit times we could not obtain 
results that agreed with the literature regarding the size of the mitotic pool or its production. 
Furthermore, a fit between our simulation model's results regarding PMN counts in peripheral 
blood with empirical data could not be achieved without speculating extensively on the nature 
of G-CSF effects on non-committed stem cells (data not provided). 

The effects of G-CSF on this compartment are modeled as an increase in the rate of cells 
entering the myeloblasts from the uncommitted stem cell pool (Eq. 1.2), increases in the rates 
of mitosis and introduction of new points of amplification. Since little data is available 
regarding the increases in amplification due to G-CSF, we first assumed that amplification 
reaches full potential at points that under normal conditions undergo an amplification of 
below a factor of 2. Additionally we assume that the transit time in all mitotic sub- 
compartments and the typical cell cycle duration are not affected by G-CSF, since as far as we 
know, there are no statistically robust data in the literature. Lord's (1989) results in this 
connection are ambivalent. 
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These modifications yield increased neutrophil production in accordance with the neutrophil 
counts reported by Price et al (figs 4,5). In addition, these increases are in good accordance 
with Price's data about neutrophil bone marrow pool sizes. 
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To the extent of our knowledge, there is no experimental data regarding the effects of G-CSF 
on the commitment of stem cells to the mitotic compartment. However, we were not able to 
reproduce the effect of G-CSF on neutrophil counts and the mitotic compartment sizes 
beyond day 5 of administration without assuming an increase (15% with the highest dose of 
G-CSF) in the rate of cells entering the myeloblast compartment. We note however, that 
alternatively, G-CSF may change the behaviour of the myeloblast compartment such that 
some of the cells there undergo self-renewal instead of moving on to the promyelocyte 
compartment. If more concrete experimental data emerges about G-CSF effects on stem-cells 
in-vivo that points otherwise the model can be modified. 

Post Mitotic Compartment 

The post mitoti c compartment doesn't inc lude m echanisms of amp ljfiraijn^^ and 
is apparently relatively insensitive to cytotoxic chemotherapy. Therefore, it is biologically 
acceptable and computationally sensible to model this compartment as a single pool of cells, 
such that the last hourly cohort of the mitotic compartment enters it, and a proportion of the 
cells in it enter the neutrophil pool every hour (Eq. 1.5). 

The effects of G-CSF on the post-mitotic compartment are shortening its transit time 
(increasing l 3 in Eq. 1.5). Price notes, that the number of cells in the post mitotic 
compartment is not significantly changed following administration of G-CSF. This 
determination is based on counts made on day 5 after G-CSF administration. We can 
therefore safely assume that any increased production of the mitotic compartment flowing 
into the post mitotic compartment is translated over the long-term to an increase in the flow of 
cells from the post mitotic compartment to the neutrophil pool. This increased flow is 
compensated by increased production in the mitotic compartment only at a later stage. 
Therefore, we set an upper limit to the number of cells in the post mitotic compartment, 
which is at the values given as steady state counts (Yl). 

In brief, the effects of G-CSF on the neutrophil lineage are modeled during the first few days 
primarily as a decrease in the counts of the post-mitotic compartment, which is then 
compensated by an increased production in the mitotic pool. This compensation sustains the 
increase in neutrophil counts in peripheral blood. Though no empirical data is available on 
this point, simulations of our model predict (fig 3.) that the number of cells in the post-mitotic 
compartment decreases substantially during the first two days of G-CSF administration, and 
then replenishes, so that on the sixth day the counts return almost to their normal levels. This 
replenishment lags behind that of Price et al report by a few hours. We can thus formulate a 
testable hypothesis, i.e., whether using the same G-CSF protocol Price et al used, there is 
indeed a nadir on day 3 of the treatment. 
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Compartment 


Day 0 (no 


Day 15 of G-CSF 


Relative increase 




G-CSF) 


treatment (xlO 9 


in compartment 




(xlO 9 cells/ 


cells/ kg. body 


size due to G- 




kg. Body 


weight) 


CSF 




weight) 






Myeloblasts 


0.140 


0.153 


1.09 


Promyelocytes 


0.582 


0.898 


1.54 


Myelocytes 


1.373 


3.564 


2.60 


Mitotic Total = 


2.10 


4.615 


2.20 


Circulating neutrophils 


0.4 


2.35 


5.88 



Table 1. Simulated kinetics after 15 days of subcutaneous administration of 300ug G-CSF/kg weight. Day 0 
values are the mean values Dancey et al (1976) use. 




Days 

Fig 3. The simulated effects of daily administration of 300[ig of G-CSF on the post mitotic 
compartment. Note the nadir on day 3 is almost entirely compensated by day 6. 



Neutrophils and G-CSF in the Circulating Blood 

The elimination of neutrophils from peripheral blood follows a Poisson distribution, 
and can therefore be described as an exponential function (Cartwright, 1964). 
Therefore the rate of cells leaving this compartment is based on half-life determinations 
available in the literature (Eq. 1.7). Since no direct cytotoxic effects of chemotherapy have 
been described for this compartment it is also modeled as a single pool of cells. 

The kinetics of G-CSF is also modeled as an exponential distribution with a half-life of 3.5 
hours (Stute, 1992) (Eq. 1.1). 



1.2 
1 



The effects of G-CSF on the kinetics of the neutrophil lineage appear not to be a linear 
function of G-CSF administration levels. Since data provided in the literature (Chatta 1994) 
only refers to two doses (30, 300(igram/kg. Body weight) we can only speculate on the effects 
of other levels of G-CSF. After trial and error analysis, we found that assuming that the 
effects of the 300-|igram dose are the maximal, at the 30-jigram its effects are about 30% of 
the maximum. The effects as a function of G-CSF level are connected piece-wise linearly. 
This way, we obtain the neutrophil levels observed clinically under both the 300 and the 30- 
jigram protocols (fig. 4). 

Discussion 

Our model successfully predicts the kinetics of the neutrophil lineage of healthy humans 
following different protocols of G-CSF. This was determined based on actual experimental 
data. 

The formulation of the model in piecewise linear terms will allow us to use this model as a 
clinical tool in three ways: Firstly, the model will determine the effectiveness of various 
protocols suggested by clinicians prior to their actual use on human patients. Secondly, the 
software will allow us to compute the optimal protocol in a given situation of neutrophil 
counts, so that the neutropenic period following chemotherapy is either shortened or 
completely avoided at a minimal cost and exposure to G-CSF. Thirdly, the model will serve 
as a constituent in a broader framework of clinical tools that will compute the most optimal 
treatment plan for chemotherapy and growth factors. 

We expect all these to help clinicians administer more rational treatment to their patients 
minimizing both suffering and medical costs. 



Amplification at 
the exit 


Mean transit time 
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Compartment 
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0.139* 
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0.558* 


Promyelocytes 
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48" 




1.4 


Myelocytes 


1 


160* 


5.84 


Post mitotic 


0 


10.96* 


* 

0.4 


Neutrophils 



Table 1. Kinetics under steady state conditions in healthy humans. 

*Dancey, 1976. 
+ Dresch, 1986. 

"Calculated based on the steady state assumption as elaborated in the main text. 
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APPENDIX C 



Optimizing cytotoxic drug delivery (administration/efficacy) for 
cancer patients 

Introduction 

Cancer is the second leading cause of mortality in the US, resulting in approximately 
550,000 deaths a year. There has been a significant overall rise in cancer cases in 
recent years, attributable to the aging of the population. Another contributing factor to 
the rise in the verifiable number of cases is the wider use _of_screening tests, such as 
mammography and elevated levels of prostate specific antigen (PSA) in the blood. 

Neither better detection nor the natural phenomenon of aging, however, can entirely 
explain the increase in new cases of tumors. Meanwhile, other cancers, like brain 
tumors and non-Hodgkin f s lymphoma, are becoming more common. Their increase 
could refl ect changes in exposures to as yet unidentified carcinog ens. Current trends 
suggest that cancer may overtake heart disease as the nation's no. 1 killer in the 
foreseeable future. As gene therapy still faces significant hurdles before it becomes an 
established therapeutic strategy, present control of cancer depends entirely on 
ch emo therapeutic methods. 

Chemotherapy is treatment withjdrugs to destroy cancer cells. There are more than 50 
drugs that are now used to delay or stop the growth of cancer. More than a dozen 
cancers that formerly were fatal are now treatable, prolonging patients' lives with 
chemotherapy. 

Treatment is performed using agents that are wi dely non-cancer-specific, killing- Cells 
th at have a high proliferation rate. Therefore, in addition to the malignant cells, most 
chemotherapeutic agents also c ause severe side-effects because of the damage 
inflicted on normal body cells. Many patients develop severe nausea and vomiting, 
become very tired, and lose their hair temporarily. Special drugs are given to alleviate 
some of these symptoms, particularly the nausea and vomiting. Chemotherapeutic 
drugs are usually given in combination with one another or in a particular sequence 
for a relatively short time. 

Chemotherapy is a problem involving many in teracti ve^ online^^proc^ses, which 
operate on different organizational levels of the biological system. It usually involves 
genomic dynamics, namely, point mutations, gene amplification or other changes on 
the genomic level, which may result in increasing virulence of the neoplasia, or in the 
emergence of drug resistance. Chemotherapy may affect many events on the cellular 
level, such as cell-cycle arrest at different checkpoints, cell transition in and out of the 
proliferation cycle, etc. Chemotherapy may also interfere with the function of entire 
organs, most notably, with bone marrow blood production. In recent years molecular 
biology and genetics has made an important step forward in documenting many of 
these processes. Yet, for assessing the contribution of specific molecular elements to 
the great variety of disease profiles, experimental biology must be provided mthLtools 
th at allo w, a JormaL and systematic ^ anajy^ ^ oLihe^intricate- interactioa between .the 
genomic, cellular and cell population processes in the host and in the disease agent. 



This system is so complex that there is no intuitive way to know how small changes in 
the drug protocol will affect prognosis. But in spite of this intricacy, attempts to 
improve chemotherapy have been carried out by "trial and error" alone, with no 
formal theory underlying the application of specific drug schedules. Such an approach 
"is apt to result in no improvement, only discouragement and little useful information 
for future planning" (Skipper, 1986). 



1. General Description 

The product is a method for identifying clinically feasible and efficacious optimal 

^— " '** i i mi mi -m in | mm] mm mMl ^ m _ ,- ltl _. 

tr eatment protocols f or anticancer drugs. We do this by using (1) optimization ~* 
methods taken from Operation Research, a nd involving heuristic search algorithms. 
These are applied for (2) testing the performance of a very large number of drug 
protocols on (3) "virtual cancer patients," in the form _of realistic computer 
simulations of detailed mathematical models of the involved dynamics. Our method 
is tffiTflrst which can identify optimal treatment protocols, without any 
constraint on the level of complexity of the mathematical models of the 
biological, clinical and pharmaceutical processes. This means that the 
description of the system can be as realistic as required. Previous optimization 
methods imposed constraints on the level of complexity of the mathematical methods, 
so that the simulated systems were always too simplistic to retrieve the clinical 
situation. 

The assumption that underlies this method is, that the effect of a certain drug on a 
certain cancer cells population depends on the overall population dynamics of the 
tumor. To describe these complicated dynamics we use a set of models thaUbllow 
up the changes in canc e^cdj^numbers, spa tial distri bution and other relevant 
phenotypic (epigenetic) factors, taking into consideration also the g enetic profile of 
the patient , in cases where this one is available. To also allow for the use of cell- 
cycle-phase-specific drugs, we follow-up each drug susceptible body cell (normal and 
cancer) according to its cell-cycle status. Because any realistic description of such 
"virtual patient" is too complex to be attacked by analytical optimization method, then 
our introduction of the heuri stic optimization methods approach is a highly sig nificant 
novelty. 

The project consists of 6 successive stages: 

1) D eveloping the appropriate mathematical mod els to describe the cancer cell 
populations with respect to the c ell-cycle as well as spatio-temp oral dynamics. This 
model can be as complex as necessary for enabling an accurate computation of the 
number of cancer cells at any given moment. 

2) Inp utting individual parameters i nto the model and, hence, bridging between the 
model and accessible clinical data regarding profiles of the tumor and of the patient. 

3) I ntroducing dru g characteristics^ (mode of action, effect and on tumor cells etc.), 
the pharmacokinetics and the pharmacodynamics into the model, using experimental 
data about the relevant drugs. 

4) Combining this mod el with additional models that con sider toxic effects of the 
drug on the host (e.g., granulocytic and thrombocytic lineages). 



5) Defining a general fitness function f or given protocol, so that this system may 
serve as prediction tool for anticancer treatment (with an option of tuning by the 
user). 

6) Using search algorithms.arid jc&^ to search through the possible 
protocols space in order to obtain the best/most suitable solution. 

The search process may be performed on several levels, starting from the most 
general (for a given type of disease), to the most detailed and patient-specific level 
(which will require relevant clinical data about the patient). In the case of a more 
generalized search, the certain region of protocols may be proposed as most suitable 
for given type/profile of cancer; in the case of specified search, this region can be 
narrowed according to additional patient data. 



2. Mathematical Model: 

2.1. Basic Tumor Layer 

The basic layer of the model incorporates a description of age distribution of cycling 
cells and number of resting (quiescent) cells. The term "age of the cell" here refers to 
chronologica Tage starting from the conven tional beginning point ofmitotic cycle. 

In what follows we give the definition of this basic model: 

- The whole cycle is divided into 4 compartments, or stages (G\> S, G 2 and M). Each 
compartment is divided into equal subcompartments, where / subcompartment in 
each stage represents cells of age i in this stage (i.e. they have spent i time-steps in 
this stage). 

- Let {G;, 5, G2, M} denote the corresponding stages of the cell cycle, and Go - the 
quiescence. 

- Let Tk denote the maximum duration of k stage in the cycle. 

- We denote as • t the time resolution of the model - discrete time steps. 

- x[ (t) is a function which represents the number of cells in stage k with age i into 
this stage, at time t to f +• t. Both time and age are measured in same units (hours). 

- Q(t) represents the number of resting cells at time t to f. 

- Trans(k,i,t) represents the probability that a cell of age i in the stage k will move to 
the next (k+1) compartment. Cells entering the new stage always start from the first 
subcompartment, i.e. from i=l. This probability may change with time. 

By definition: V& : Trans(k,T k ) = 1 That means that the cell cannot remain more than 
Tic in k compartment. 

- Let To and Tr denote the onset of Gl and the age at which the cell passes through R- 
point (restriction point) in Gl, respectively. The later is associated with commitment 
to complete the mitotic cycle. 

- The total number of proliferating cells P(t) can be calculated as follows: 




- In every time interval, quiescent cells may return to the proliferation pool. 
Alternatively proliferating cells may change their state to become quiescent if they 
are in Gl stage and at age i, where Tr* i>0 . To describe this process we introduce the 
function G;. o(i,t) which describes the number of Gl cells in age i which become 
quiescent , during time interval [t, r+ • t]. This function may receive negative 
values, accounting for cells that return from resting to proliferation. 

As we assume that the exit to quiescence can occur only prior to the R-checkpoint 
(even in cancer cells), and that a resting cell which returns to proliferation enters the 
cycle at To, it can be stated: 

\/i>T R ,Vt:G^ 0 (i,t) = 0 

It must be noted, that this function is not dependent on i and t solely. Its value is 
determined according to current cell distribution and all the general parameters that 
characterize the described cells group. The same should be said about the values of 
Trans vector that can change during the history of given population. 

- The model traces the development of described group of cancer cells using given 
parameters, by calculating the number of cells in each and every subcompartment 
according to the following stepwise equations: 

x M (/-!)■[!- Transik, i - 1, t - 1)] , 1 < i < T k 
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for k=G2,M, k-1 returns the previous stage (e.g. G2-I-S). 
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- These equations make it possible to calculate the number of cells in each 
subcompartment at every time interval [t* t] starting from initial distribution (e.g. at 
time t-0). 

Remarks: 

1) As stated above, cells ages are measured in absolute time units, and therefore refer 
to the chronological age of the cell (as opposed to biological age, which refers to 
maturity of the cell and its units are relative to maturation rate that is different from 
cell to cell). As a consequence, in this model no cell can remain in the same age 
subcompartment after every time step. On the other hand, a fraction of the cells that 
leaves any subcompartment may be transferred to the first subcompartment of the 



next stage, according to probability vector Trans(k,i,t). This fact enables us to account 
for variability of cycle lengths while retaining deterministic approach. 

2) The behavior of the cell population in this model is wholly controlled by two 
components: Trans vector, and GltoGO function. These two functions determine 
uniquely the outcome of every single time step, and, consequently the result over the 
long periods. We name them here "control functions". As we have mentioned, the 
values of these functions may be dependent not only on time and age of cells, but also 
on the current population state (or, generally, on the whole history of the population) 
as well as on the environment associated with given cells group. On the other hand, 
those parameters are similar for all the cells in the group. This implies that the model 
presented here is suitable for describing highly homogenous group of cells. As it will 
be claimed later, control functions values are determined by the parameters that 
describe the local conditions for given population. Therefore, we expect the basic 
layer of the model to give realistic description of uniform group of cancer cells for 
which environmental conditions and relevant biological properties are defined, in a 
way that will allow the construction of the control functions for this group. 

2.2. The General Tumor Model 

In the general approach the whole model is viewed as constructed from similar 
components, each of them derived from the basic structure described in the previous 
paragraph. Each component represents cells that are situated in the same 
environmental conditions and, therefore behaves similarly. We call them groups. The 
whole tumor is modeled as union of heterogeneous groups of cells, where the 
development of every group can be accurately predicted (when local conditions are 
known). 

This general model simulates progress of tumor in time by discrete steps. At each 
step the number of cells in each subcompartment of each group is calculated 
according to the previous state, parameters of tumor, drug concentration, etc. The 
parameters of the tumor must include all the information that is relevant to the 
prognosis. Some of these parameters are defined locally, e.g., those relating to the 
tumor 's geometry. For this reason the representation of the spatial structure will be 
included. 

The cells will be able to pass between the groups during the development of the 
tumor. This allows the representation of the changes in the local conditions, during 
the tumor evolution (e.g. forming of necrotic core, improvement in "living 
conditions" in vascularized regions, etc.). In addition, all the parameters of the tumor 
may change in accordance with the dynamics of the cancer. 

The calculation of the tumor development over time will be done by stepwise 
execution of the described simulation and can be used to predict the outcome of the 
treatment or in fitness function for search algorithms. 

2.3, From General to Individual Tumor Model 

When the general theoretical description of the model is accomplished, the model is 
fitted to represent the actual tumors. We render it patient-specific by adjusting all the 
parameters that determine the behavior of the modeled tumor to those of the real 
cancer in the patient's body. In order to accomplish this task we will establish the 



connections between mathematical parameters (most of them will have direct 
biological implication) and every kind of data that is practically obtainable in the 
clinic. This individualization procedure includes three aspects: 

1) individual parameters of tumor dynamics 

2) individual parameters of specific drug pharmacokinetics for a specific patient. 

3) individual parameters of the dynamics of dose-limiting normal host tissues. 
These connections may be defined through research on statistic correlation between 

different parameters (including genotype-phenotype correlation), or using advanced 
biochemical research (which may establish quantitative relations). 

In the general model important dynamic parameters are estimated from experimental 
studies that are conducted in certain patient populations. Any of these parameters, 
when available on the per patient basis, can be individualized, while others left as 
population-based means. This approach allows continuous increase in the degree of 
individualization of the treatment protocols with progress in the technology. 

Defined so, the model will be able to give realistic predictions for treatment outcome 
either for specific patient or for a broad range profiles of patients and diseases. This 
tool can serve to perform the prognosis of an untreated cancer patient, or as a basis for 
treatment modeling. 

3. Introducing Pharmacology 

In order to simulate cancer treatment we add to this model the pharmacologic 
component. We model pharmacokinetics as well as pharmacodynamics for specific 
anticancer drugs. We begin by representing cell-cycle specific drugs, however the 
model implies no restriction on the type of drugs to be employed. 

We model the distribution of the drug in and around the tumor as well as in the 
blood (the drug kinetics). For this purpose, we use the suitable model, defining it 
precisely for every certain type of the drug. The concentrations of drug in the body are 
calculated at every time step in accordance with the drug administration specified by 
the protocol. 

The dynamics of the drug are represented through the direct influence of the drug 
on tumor cells. The effects on the proliferating cells are mostly blocking the cycle in 
different stages (which can be modeled as cell arrest) and cell death (immediate or 
after being in the block). Cell-cycle specific drugs have no direct influence on 
quiescent cells, but can affect them indirectly by killing proliferative cells and 
therefor changing local conditions. Where additional types of drugs added to the 
model, their effect on any kind of cells may be too modeled as killing certain fraction 
of cells (which is dose-dependent) or changing the behavior of the cells. 

Additional phenomena that may come out to be significant in drug kinetics and 
dynamics (e.g. rate of absorption by the cells, development of tumour resistance to 
specific drug, etc.) can be introduced into the model to make it as realistic as needed. 

The description of the drug in the model is done in terms of quantitative functions, 
which enable to calculate the drug amounts at certain locations and the tumor 
response to it at every time step. In the general case, these functions include 
parameters that depend on the specific data (drug type, body parameters, 
characteristic of the tumor, etc.) and can be determined in given situation. The relation 
between clinical data and these parameters can be established in the ways similar to 
those described for the cancer model. 

The combination of cancer model with the drug model described above makes it 
possible to predict the outcome of the treatment, given the relevant parameters for the 



drug, the cancer and the patient. Again, the prognosis may be made for specified cases 
as well as for broad profiles of patients or disease. This simulation also serves to build 
the fitness function used in search algorithms in optimization system. 

4. Combining with minimizing host toxicity. 

Although being accurate predictive tool, the model that represents chemotherapy of 
tumor alone cannot be used in optimization, for it posts no constraint on the choosing 
of the treatment. Actually, this model implies using as much drug as it is possible 
until the final elimination of the tumor; while in the live system the toxicity of the 
drug is the most important constraint limiting the treatment. In most cases of 
anticancer chemotherapy the dose-limiting toxicity is bone marrow suppression, the 
two most sensitive bone marrow lineages being granulopoiesis and thrombopoiesis. 
Accordingly, those two were chosen as an example and are modeled separately and in 
the similar way to predict the negative effect of the chemotherapy on them. These 
models reconstruct the damage caused by the chemotherapy to the bone marrow cells 
and the recovery of these lineages (also treated by specific growth factors). 

Thus, the whole system is able to predict the result of chemotherapy treatment for 
the tumor as well as for bone marrow cells, allowing the use of the protocols that 
combine anticancer drugs and growth factors for healthy cells. 

Chemotherapy toxicity to any other normal host cell population can be similarly 
taken into account, if it is proven to be relevant for dose and schedule optimization. 



5. Defining the Fitness Function 

Fitness function is an important tool in operation research. In the method developed 
by us it allows the comparison between a number of different protocols each one of 
them scoring differently with respect to various objectives that can be set by the 
developers or by the users and identifying the protocol for which we predict the best 
weighted score. Fitness function calculates for any given protocol its relative efficacy, 
thus enabling a definite decision of the best protocol from a given set of protocols. 

In different cases, different objectives can be formulated. There are 3 main settings 
in which such a model can be used: 
One) clinical practice- where objectives can change depending on type of disease, 

condition of the patient, purpose of treatment, etc; 
Two) pharmaceutical company- where objectives can be aimed at finding the 

therapeutic window and an optimal schedule; 
Three) scientific setting- where different and sometimes esoteric objectives 

can be aimed at. 

For any such case, certain fitness function can be formulated, reflecting all given 
requirements. Thus, in any particular case one can compare between different 
protocols and obtain the most suitable to his/her special purposes and needs. 

Using the fitness function we can: a) estimate the efficacy of a given protocol, b) 
search for the solution of an optimization problem, that is predict which protocol will 
be best of many potential protocols considered for curing/relieving the patient. 



6. Solving the Optimization Problem 



The optimization problem is stated using the described models: to find the protocol 
for drug administration (with option to growth factor administration) which 
maximizes the given fitness function. 

As explained above, we define the fitness function according to the user 
requirements. For example, the goal of the treatment may be defined as minimizing 
the number of cancer cells at the end of the treatment, minimizing the damage to the 
BM cells throughout the treatment or at its end, and curing the patient (where cure is 
defined precisely) as quickly as possible. Note that the fitness function may also 
include goals such as maximizing life expectancy, minimizing cost of treatment, etc. 

Generally, the aim of optimization is to find the best protocol, i.e. the protocol that 
generates the best value of fitness. Customarily, this is achieved by mathematical 
analysis. However, mathematical analysis is restricted to relatively simple models, 
whereas our models are very complex and, therefore, cannot be solved analytically. 
On the other hand, the practical purpose of the treatment is not to find the best 
possible protocol (i.e., the global optimum) but "only" one that will suit the user's 
objectives, even if its fitness is not absolutely the best (i.e., the local optimum). For 
this reason we can be satisfied with the solution that can be shown to promise that the 
patient will be cured. 

Hence, the optimization problem may be reformulated as follows: for given starting 
conditions, find the treatment protocol that will fulfill the user's requirements (e.g. 
curing a patient according to given definitions of cure) and subjected to given 
limitations (e.g. treatment duration, drug amounts, etc.). To this end it is unnecessary 
to find the global solution. It is enough to perform search, using search algorithms, in 
certain regions of the protocols' space, and find the local maxima of the fitness 
function. After determining these locally best protocols, we can verify that they serve 
one's objectives and check them numerically for stability. 

Such a strategy can be used for identifying patient-specific treatment, as well as in 
the general case, where only the profile of the disease and the drug specified. If more 
patient-specific data are supplied, the solution will be tailored more specifically. On 
the other hand, the optimization program could propose general recommendations for 
the protocol types for certain kinds of tumor. 
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Optimizing drug regimens in cancer chemotherapy by an 
efficacy-toxicity mathematical model. 

Iliadis A, Barbolosi D. 

Faculty of Pharmacy, University of 
Marseilles, 27, boulevard Jean Moulin, Marseille Cedex 5 13385 
France. 

In cancer chemotherapy, it is important to design treatment strategies 
that ensure a desired rate of tumor cell kill without unacceptable 
toxicity. To optimize treatment, we used a mathematical model 
describing the pharmacokinetics of anticancer drugs, antitumor 
efficacy, and drug toxicity. This model was associated with 
constraints on the allowed plasma concentrations, drug exposure 
and leukopenia. Given a schedule of drug administrations the 
mathematical model optimized the drug doses that can minimize the 
tumor burden while limiting toxicity at the level of the white blood 
cells The main result is that the optimal drug administration is an 
initial high-dose chemotherapy up to saturation of constraints 
associated with normal cell toxicity and a maintenance continuous 
infusion at a moderate rate. Data related to etoposide investigations 
were used in a feasibility study. Simulations with the optimized 
protocol showed better performances than usual clinical protocols 
Mode -based optimal drug doses provide for greater cytoreduction 
while limiting the risk of unacceptable toxicity. Copyriqht 2000 
Academic Press. 
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